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SECTION  I 
INTRODUCTION 


A. 


THE  TASK 


Source  identification  of  seismic  events  has  been  studied  for 
two  decades.  As  a rule,  these  research  efforts  have  evaluated  discrimination 
techniques  applied  to  events  from  a particular  region  as  recorded  at  a partic- 
ular station.  The  drawback  to  this  type  of  study  is  that  it  fails  to  consider  the 
effects  of  source,  path,  and  receiver  variations  on  the  discriminant.  Thus, 
identification  techniques  such  as  first  motion  or  complexity,  which  initially 
appeared  to  be  quite  promising,  were  later  dropped  or  downgraded  in  im- 
portance when  considered  in  light  of  larger  data  bases.  Even  a technique 
such  as  Mg-  m^,  which  has  historically  been  shown  to  be  highly  effective  in 
discriminating  between  earthquakes  and  nuclear  explosions,  must  be  applied 
With  care.  For  example,  Peppin  and  McEvilly  (1974)  used  Pn  versus  12  sec- 


ond LR  (analogous  to  M - m ) to  discriminate  between  Nevada  earthquakes 

8 D 


and  NTS  explosions.  Their  plots  of  data  with  ML  less  than  5.  0 showed  good 
separation  between  the  two  populations.  However,  inclusion  in  their  plots  of 
explosion  Pn  versus  LR  points  for  greater  than  5.  0 shows  that  these  ex- 
plosion points  fall  in  the  earthquake  population.  Thus,  the  discriminant  fails 
for  this  region  when  considering  large  magnitude  events. 


Since  no  one  discrimination  technique  can  be  expected  to  cor- 
rectly identify  all  events  from  all  regions  as  recorded  at  all  stations,  it  is 
necessary  at  this  time  to  attempt  to  identify  events  using  all  available  dis- 
crimination information  in  a multivariate  mode  of  operation.  The  utility  of 
this  type  of  approach  is  demonstrated  by  Anglin  (1971),  who  shows  that  a 
complete  separation  of  a suite  of  Eurasian  earthquakes  and  nuclear  explosions 
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could  be  achieved  by  plotting  complexity  discriminant  values  against  third 
moment  of  frequency  discriminant  values.  Used  individually,  neither  of 
these  discriminants  completely  separated  the  two  populations. 


The  specific  goals  of  this  report  are  therefore: 

• To  define  a set  of  identification  criteria  which  can  be  applied 
to  network  data  to  optimally  use  the  features  of  each  discrim- 
ination method. 

• To  determine  the  physical  source  of  each  discrimination  meth- 
od in  order  to  allow  the  application  of  the  method  to  events  of 
similar  physical  characteristics. 

• To  assemble  a package  of  discrimination  techniques  and  test 
this  package  in  order  to  define  how  the  techniques  perform, 
both  individually  and  in  a multivariate  mode. 

• To  recommend  further  package  refinements  based  on  any  de- 
ficiencies or  operational  difficulties  discovered  during  the 
course  of  this  study. 

i 
1 

The  most  significant  progress  toward  these  goals  last  year 
was  the  development  of  an  integrated  event  identification  system.  The  purpose 
of  this  system  was  to  automate  as  fully  as  possible  the  routine  stages  of  the 
data  processing  so  that  the  analyst  can  devote  as  much  time  as  possible  to 
those  areas  which  resist  automation.  In  brief,  the  system  consists  of  two 
stages  - measurement  and  analysis.  The  measurement  stage,  which  is  doc- 
umented by  Schmidt  and  Wilson  (1978),  edits  the  data  from  the  event  tape, 
removes  the  instrument  response  to  produce  displacement  data,  and  measures 
the  parameters  to  be  used  in  the  discrimination  analysis.  This  stage  of  the 
system  is  performed  on  the  IBM  360/44  computer.  The  edited  waveforms 
and  measured  parameters  are  written  on  tape  for  transfer  to  the  PDP-15/50 
computer.  At  the  PDP-15/50,  the  analyst  is  able  to  visually  inspect  the  data 
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in  a quality  check  mode,  override  previously  made  measurements  as  desired 
compute  unbiased  network  magnitudes  using  Ringdal's  technique  (Ringdal, 
1975),  and  perform  multivariate  discrimination.  All  of  this  will  be  more 
fully  discussed  in  later  sections  of  this  report. 


It  is  important  to  note  that  the  work  carried  out  during  the  last 
year  had  as  its  principal  goal  the  creation  of  this  system,  since  the  long-term 
benefits  of  such  a system  when  fully  operational  will  far  outweigh  the  short- 
term benefits  of  measuring  a number  of  discriminant  values  on  any  particular 
data  set.  The  discrimination  methods  embodied  in  the  current  system  are  not 
intended  to  represent  all  possible  methods.  The  intent  in  the  work  performed 
during  the  last  year  was  to  build  in  enough  recognized  discrimination  techniques 
to  permit  measurement  of  a data  set  adequate  for  the  multivariate  discrimina- 
tion effort.  By  modularizing  the  system,  it  will  always  be  possible  to  add 
more  discriminant- measuring  methods  at  a later  date. 


An  important  point  to  note  is  that  past  research  has  indicated 
that  discriminant  measurements  are  strongly  influenced  by  the  tectonic  re- 
gion of  the  source.  Past  efforts  to  minimize  the  influence  of  the  source  tec- 
tonic region  have  generally  taken  the  form  of  separating  the  events  by  geo- 
graphic region.  For  example,  the  Mg-  discriminant  using  SRO  data  was 
presented  in  this  fashion  by  Strauss  and  Weltman  (1977).  However,  this  ap- 
proach is  not  optimal,  since  when  studying  a large  area  such  as  Eurasia  many 
regions,  and  hence  many  test  events,  must  be  studied.  In  an  attempt  to  sim- 
plify thi6  problem,  the  events  are  classed  as  being  located  in  tectonically 
active,  inactive,  or  rift  regions.  Since  the  tectonic  region  in  which  the  sta- 
tion is  located  can  also  influence  discriminant  measurements,  the  stations 
are  also  classified  as  being  located  in  tectonically  active,  inactive,  or  rift 
regions.  This  results  in  a total  of  only  nine  station- event  region  tectonic 
classifications.  Figure  1-1  illustrates  the  world-wide  distribution  of  the 
three  tectonic  regions. 
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THE  AREA  OF  INTEREST  DATA  BASE 
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The  data  base  to  be  used  in  this  event  identification  study  is 
the  suite  of  'Area  of  Interest'  (AI)  events  supplied  to  all  researchers  by  Tele- 
dyne Geotech  Incorporated  of  Alexandria,  Virginia.  The  data  supplied  for 
each  event  consists  of  edits  of  the  data  recorded  at  the  stations  listed  in 
Table  1-1  and  shown  in  Figure  1-2.  (An  asterisk  in  Table  1-1  denotes  that 
the  station  was  used  in  this  research  effort.  ) These  edits  nominally  consist 
of  180  seconds  of  short-period  P,  600  seconds  of  short-period  Lg,  240  sec- 
onds each  of  long-period  P and  S,  and  1200  seconds  of  long-period  Rayleigh, 
where  the  edit  windows  are  roughly  centered  about  the  expected  arrival  time 
of  the  appropriate  waveform.  The  only  change  made  in  the  data  when  pre- 
paring the  AI  tapes  was  to  convert  the  data  from  digital  counts  to  millimicrons. 
The  data  from  all  stations  except  RKON  and  HNME  are  in  the  vertical,  north, 
east  configuration.  The  data  from  RKON  and  HNME  were  recorded  in  a ver- 
tical, transverse,  radial  configuration  with  respect  to  the  Nevada  Test  Site. 


As  of  30  September  1978,  the  requisite  parameters  (date,  or- 
igin time,  and  epicenter  coordinates)  for  87  AI  events  had  been  supplied  to 
the  authors.  Since  the  event  identification  system  as  presently  constituted  is 
designed  for  the  shallow  earthquake/explosion  discrimination  problem,  it 
was  necessary  to  eliminate  from  the  data  base  all  events  with  depth  unknown 
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or  greater  than  60  km  by  searching  the  National  Earthquake  Information  Ser- 
vice bulletin  for  depth  information.  This  depth  information  requirement  is 
only  a temporary  measure  to  permit  testing  of  the  system  as  presently  de- 
signed. Necessary  future  work  on  the  system  will  include  placing  a depth- 
determination  algorithm  such  as  that  designed  by  Page  (1976)  to  provide  a 
screening  procedure  (Basham  and  Anglin,  1973)  which  will  eliminate  most  of 
those  events  with  focal  depths  greater  than  the  maximum  feasible  drilling 
depth  after  making  allowance  for  uncertainty  in  depth  determination. 


TABLE  1-1 

STATION  LOCATIONS 


Site  No. 


Station 

Latitude 

Longitude 

ANMO* 

34.  94 

-106.  46 

ANTO 

40.  00 

33.  00 

BOCO 

5.  00 

-74. 00 

CHTO* 

18.  80 

99.  00 

NORS* 

60.  84 

10.  89 

GUMO* 

13.  59 

144.  87 

MAIO* 

36.  30 

59.49 

LASA* 

46.  69 

-106.  22 

NWAO* 

-32.  93 

117.  24 

KSRS* 

37.  45 

127.  92 

SHIO 

25.  57 

91.  88 

TATO* 

25.  00 

121.  50 

SNZO* 

-41.  31 

174.  70 

ILPA* 

35.  70 

50.  61 

ALPA* 

65.  03 

-147.  20 

CTAO  (ASRO)* 

-20.  09 

146.  25 

ZOBO* 

-16.  27 

-68.  12 

KAAO* 

34.  54 

69.  04 

MA  JO* 

36.  54 

138.  21 

ATAK* 

52.  88 

173.  17 

BFAK* 

64.  77 

-146. 89 

CTAO  (HGLP)* 

-20. 09 

146. 25 

CHGO 

18.  79 

98.  98 

TNAK* 

62.  91 

-156. 02 

TLOO 

39.  86 

-4.  01 

EIAO 

29.  55 

34.  95 

KONO 

59.  65 

9.  60 

OGDO 

41. 07 

-74.  62 

KIPO 

21. 42 

-158.  02 

ALQO 

34.  94 

-106.  46 

ZLPO 

-16.  27 

-68.  12 

MATO 

36.  54 

138.  21 

• HNME* 

46.  16 

-67. 99 

RKON* 

50.  84 

-93.  67 

t A s Active,  I = Inactive,  R 
* Stations  used  in  this  study. 
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Since  it  is  desirable  to  have  as  clean  a data  set  as  possible  be- 
fore making  the  trial  run  of  the  system,  the  plots  of  the  AI  data  supplied  by 
Teledyne  Geotech  Incorporated  were  visually  reviewed  for  mixed  events  and 
station  recording  malfunctions  such  as  spikes.  Data  showing  either  of  these 
problems  were  dropped  from  further  consideration,  since  any  discrimination 
parameters  measured  on  these  data  would  not  be  representative  of  the  event 
under  consideration. 

The  above  considerations  coupled  with  the  shortage  of  time 
available,  reduced  the  initial  data  base  to  26  events.  Upon  the  request  of 
the  contract  monitor,  nine  more  events  were  added  to  the  processing  data 
base,  bringing  the  total  to  35  events.  The  origin  times  and  coordinates  of 
these  events  are  listed  in  Table  1-2,  where  the  event  numbers  are  those  as- 
signed by  Teledyne  Geotech  Incorporated.  The  approximate  locations  of  all 
87  events  are  shown  in  Figure  1-3,  where  the  events  used  in  this  study  are 
indicated  by  their  event  numbers. 

Several  problems  were  encountered  with  the  AI  data  when  pro- 
cessing was  initiated.  First,  it  was  found  that  the  waveform  identification 
scheme  was  not  always  followed  when  the  AI  data  tapes  were  created.  Ac- 
cording to  a memorandum  issued  by  Teledyne  Geotech  Incorporated,  each 
waveform  was  to  have  a unique  seismogram  number  identifying  the  event, 
recording  station,  type  of  data  (short-period  or  long-period),  and  propaga- 
tion mode  (P,  S,  Lg,  or  Rayleigh).  However,  it  was  found  that,  for  more 
recently  created  AI  tapes,  it  was  not  uncommon  for  only  one  seismogram 
(identified  as  long-period  P)  containing  long-period  P,  S,  and  Rayleigh  to  be 
on  the  tape.  Thus,  one  cannot  with  confidence  create  a program  to  search 
on  the  basis  of  seismogram  number  for  one  propagation  mode  but  must  either 
check  header  dumps  of  the  tapes  or  have  the  program  check  each  waveform 
for  the  desired  time  window  based  on  start  time  and  length  of  the  edit  and 
expected  arrival  time  of  the  desired  propagation  mode. 
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Event 

Number 


TABLE  1-2 

EVENT  PARAMETERS 


Date 

(Mo/Da/Yr) 

Origin  Time 
(Hr:Min:Sec) 

Latitude 

(°N) 

Longitude 

(°E) 

07/26/77 

16:59:59.  9 

69.4 

90.  4 

11/01/77 

03:54:24. 0 

55.  3 

130.  8 

11/04/77 

23:54:51.  5 

30.  7 

81. 3 

07/30/77 

01:56:59.  9 

49.  7 

78.  2 

08/17/77 

20:00:00.  7 

50.  9 

111.  0 

08/17/77 

04:26:59.  8 

49.  8 

78.  2 

08/20/77 

22:00:00.  6 

64.  1 

99.  8 

09/01/77 

03:00:00.  0 

73.  3 

54.  3 

09/05/77 

03:03:00.  2 

50.  1 

78.  9 

09/10/77 

16:00:00.  5 

57.  2 

106.  8 

09/30/77 

07:00:00.  0 

48.  0 

48.  0 

11/18/77 

05:20:10 

33.  0 

89.  0 

11/18/77 

15:10:10 

28.  0 

90.  0 

11/18/77 

17:23:25 

33.  0 

89.  0 

10/16/77 

20:03:35 

48.4 

152.  9 

11/18/77 

21:55:37. 0 

60. 1 

143.  2 

10/16/77 

15:02:49 

36.  9 

71.  5 

11/18/77 

23:12:49 

33.  0 

89.  0 

10/16/77 

21:05:35 

49.  7 

155.  1 

10/19/77 

05:02:00 

36.  3 

71.  3 

10/19/77 

21:20:37 

49.  5 

155.  4 

10/20/77 

08:18:04 

56.  3 

164.  1 

10/29/77 

03:06:59.  7 

49.  8 

78.  0 

10/26/77 

05:38:52.  2 

49.0 

155.  8 

10/26/77 

07:11:31.  3 

46.  4 

153.  5 

10/28/77 

21:15:11.  5 

39.  8 

71.  9 

10/29/77 

06:26:42.  5 

41.  0 

63.  7 

10/29/77 

10:33:59.  4 

47.  3 

153.  1 

10/20/77 

23:40:35 

33.  1 

88.  1 

11/22/77 

11:33:45 

43.  0 

89.  0 

11/26/77 

22:46:46 

37.  0 

115.  0 

11/30/77 

04:06:59 

49.  9 

78.  8 

12/02/77 

12:57:10 

52.  9 

159.  7 

12/06/77 

10:52:53 

41. 4 

69.  7 

12/07/77 

16:19:33 

35.  6 

94.  5 

Class 
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LOCATIONS  OF  EVENTS 


A second  problem  was  that  in  no  case  was  there  a specific 
short-period  Lg  edit  available  from  any  of  the  SRO  or  ASRO  stations.  Al- 
though the  short-period  automatic  detector  installed  at  these  stations  must 
be  held  responsible  in  most  cases  for  this  lack  of  Lg  data,  one  would  expect 
Lg  edits  for  some  of  the  SRO-ASRO  recorded  events,  since,  as  shown  by 
Strauss  and  Weltman  (1977),  at  some  stations  the  automatic  detector  after 
being  triggered  will  often  stay  on  for  long  periods  of  time.  Thus,  one  would 
expect  there  to  be  at  least  an  occasional  Lg  edit  for  these  stations  on  the  AI 
tapes.  This  unfortunate  lack  of  Lg  edits  also  occurred  for  data  recorded  at 
the  stations  which  continuously  record  short-period  data.  It  was  not  uncom- 
mon to  have  a short -period  P wave  edit  but  no  Lg  edit  where  one  would  ex- 
pect the  data  to  have  been  available. 


A surprising  shortcoming  of  the  AI  data  base  is  the  absence  of 
short-period  S wave  edits.  Occasionally,  an  S wave  can  be  found  in  the  P 
wave  edit  when  the  station  is  near  the  event  epicenter.  However,  no  specific 
S wave  edits  were  made  during  the  completion  of  this  data  base.  This  lack 
of  short-period  S wave  data,  coupled  with  the  scarcity  of  Lg  data,  seriously 
affects  the  quality  of  the  event  identification  effort,  since  it  minimizes  the 
utility  of  the  compressional-to- shear  ratio  type  of  discriminants  which  have 
been  proven  quite  useful  in  the  past  (Booker  and  Mitronovas,  1964). 


REPORT  ORGANIZATION 


The  following  sections  of  this  report  will  discuss  in  detail  the 
work  carried  out  during  the  last  year  toward  achieving  the  goals  of  the  event 
identification  task.  Section  II  of  this  report  presents  a discussion  of  the 
theoretical  background  of  event  identification,  followed  by  a description  of 
the  particular  event  identification  methods  used  in  this  study.  Included  in 
the  discussion  of  each  method  is  the  historical  background  of  that  method  and 
an  attempt  to  relate  the  method  to  the  previously  described  theory.  Section 


Ill  presents  an  in-depth  examination  of  the  signal  identification  parameter 
measurements  program  TISSPROG.  Particular  emphasis  is  placed  on  the 
flow  of  data  through  the  sections  of  the  program  which  measure  the  event 
identification  parameters,  since  it  is  this  portion  of  the  program  which  has 
been  developed  this  last  year.  Section  IV  presents  a discussion  of  the  analy- 
sis stage  methodology  of  this  research  effort.  In  it  is  described  the  manner 
in  which  the  network  discriminant  values  are  computed  from  the  individual 
station  discriminant  values  for  each  event.  Following  this  is  a description 
of  the  method  of  multivariate  discrimination.  Section  V presents  the  actual 
analysis  of  the  35  AI  events  processed  by  this  event  identification  system. 
Section  VI  presents  the  conclusions  reached  during  the  course  of  this  research 
effort  and  recommendations  for  possible  improvements  in  the  methodology. 
Section  VII  lists  the  references  cited  in  this  report.  Finally,  Appendix  A 
discusses  the  narrowband  filter  used  in  the  measurement  of  several  different 
event  identification  parameters. 


SECTION  II 

THEORETICAL  MOTIVATION 

The  purpose  of  this  section  is  to  provide  a theoretical  r *iva- 
tion  for  the  discriminants  used  in  this  study.  Much  effort  has  been  made  in 
recent  years  to  understand  the  earthquake  and  underground  explosion  source 
mechanisms,  and  to  determine  which  features  of  the  respective  radiation 
fields  may  be  used  to  characterize  each  type  of  seismic  event.  As  a result  of 
this  research,  several  discriminants  have  evolved  which  prove  to  be  quite 
useful  in  classifying  seismic  events  either  as  earthquakes  or  as  underground 
explosions. 

A.  THEORETICAL  EARTHQUAKE  AND  UNDERGROUND  EXPLOSION 
SOURCE  MODELS 

The  underground  explosion  source  has  been  described  by  Sharpe 
(1942)  and  others  (Bishop,  1963;  Toksoz  et  al.  , 1964;  Carpenter,  1967;Lieber 
mann  and  Pomeroy,  1969).  A shock  wave  proceeds  radially  outward  from  the 
detonation  point.  In  the  region  where  this  shock  wave  propagates,  the  peak 
stress  levels  greatly  exceed  the  material  strength  of  the  surrounding  medium, 
resulting  in  rapid  dissipation  of  energy  in  the  form  of  non-linear  deformation 
(i.  e.  heating  and  fracture).  At  some  distance  from  the  detonation  point,  the 
stress  levels  associated  with  the  shock  front  fall  below  the  elastic  limit  of  the 
medium,  and  reversible,  infinitesimal  strains  prevail.  At  this  point,  energy 
is  dissipated  in  the  form  of  seismic  radiation. 

The  transition  boundary  between  the  regions  of  non-linear  and 
linear  deformation  is  not  well  defined.  However,  it  has  been  found  useful  to 
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represent  this  boundary  schematically  as  a spherical  surface,  termed  the 
'equivalent  cavity'  (Sharpe,  1942).  The  radius  of  the  equivalent  cavity  depends 
upon  the  size  of  the  explosive  charge  and  the  material  properties  of  the  medium 
in  which  the  charge  is  detonated  (Bishop,  1963).  Using  far-field  observations, 
it  is  possible  to  recover  the  pressure  history  only  as  it  exists  at  this  boundary 
(Toksoz  et  al.  , 1964;  Liebermann  and  Pomeroy,  1969). 

The  earthquake  source  mechanism  has  been  treated  rather  ex- 
tensively in  the  seismic  literature  (substantial  lists  of  references  may  be 
found  on  this  subject  in  Brune,  1970;  Hanks  and  Wyss,  1972;  and  Turnbull, 

1976).  One  of  the  earliest  descriptions  of  the  earthquake  source  is  the  Elastic 
Theory  of  Fracture  due  to  Reid  (191 1),  based  upon  observations  made  after 
the  San  Francisco  earthquake  of  1906.  According  to  the  Elastic  Rebound  Theory, 
an  earthquake  occurs  as  the  result  of  the  accumulation  of  elastic  stresses  with- 
in a given  region  of  the  earth.  These  stresses  slowly  increase  with  time  until 
they  exceed  the  strength  of  the  medium,  at  which  point  the  accumulated  stresses 
are  relieved  by  fracture. 

The  rupture  mechanism  involved  in  the  earthquake  source  was 
originally  thought  to  be  simple  Coulomb  fracture.  However,  consideration  of 
shear  strengths  in  the  crust  and  upper  mantle,  and  the  rate  of  increase  of 
overburden  pressures  with  depth,  indicates  that  Coulomb  fracture  is  not  a 
major  contributing  factor  to  the  earthquake  rupture  mechanism  at  depths  greater 
than  25  kilometers  (Stacey,  1969). 

■ 

_ 1 

Recent  attempts  to  explain  the  physics  of  the  rupture  phenomenon 

have  resulted  in  several  fairly  successful  earthquake  models.  Steketee  (1958) 
demonstrated  the  applicability  of  elastic  dislocation  theory  to  the  earthquake 
source  problem.  In  this  case,  rupture  formation  is  initiated  by  the  sudden 
application  of  stresses  in  an  unstrained  medium,  equivalent  to  the  instantaneous 
removal  of  stresses  across  a fault  surface. 
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Dislocation  models  have  been  used  successfully  in  several 
studies  of  the  far-field  radiation  and  static  deformation  caused  by  earthquakes 
(Chinnery,  I960,  1961;  Knopoff  and  Gilbert,  I960;  Burridge  and  Knopoff,  1964; 
Haskell,  1964;  Berkhemer  and  Jacob,  1968).  Additional  studies  (Aki,  1968; 
Haskell,  1969)  have  used  dislocation  models  to  explain  near-source  displace- 
ments. A successful  earthquake  dislocation  model  has  been  developed  by 
Brune  (1970).  This  model  is  derived  by  considering  the  effective  stress  avail- 
able to  accelerate  opposing  sides  of  the  fault.  Brune's  model  provides  near- 
and  far-field  displacement-time  functions  and  spectra  for  shear  waves  which 
agree  reasonably  well  with  observed  earthquake  spectra.  In  addition,  Brune's 
model  also  provides  a possible  means  for  estimating  the  effective  stress,  stress 
drop,  and  source  dimensions  by  comparison  of  observed  seismic  spectra  to 
theoretical  spectra  obtained  from  the  model.  This  technique  has  been  used  to 
estimate  source  parameters  for  several  large  earthquakes  (Hanks  and  Wyss, 
1972). 

Attempts  to  explain  the  dynamics  of  the  earthquake  rupture 
mechanism  have  also  led  to  the  stress  relaxation  model  of  Archambeau  (1968), 
where  the  rupture  process  is  modeled  in  terms  of  a volume  relaxation  phe- 
nomenon. In  this  theory,  the  dynamical  properties  of  the  source  are  described 
in  terms  of  the  relaxation  or  readjustment  of  pre-existing  stress  fields.  The 
initial  shear  stress  field  is  taken  to  be  the  primary  cause  of  rupture.  Seismic 
effects  are  the  result  of  the  relaxation  of  this  prestress  shear  field  in  the 
medium  surrounding  a region  whose  physical  properties  (rigidity  in  particular) 
have  changed  suddenly. 

The  tectonic  source  is  viewed  in  terms  of  the  release  of  potential 
strain  energy  from  within  a non-elastic  zone  in  the  medium.  Within  this  non- 
linear rupture  zone,  stress  release  is  achieved  by  flow,  fracture,  or  phase 
change,  so  that  strain  energy  within  the  volume  is  dissipated  in  the  work  of 
non-linear  deformation. 
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This  dynamic  relaxation  model  requires  that  the  stress  field 
surrounding  the  rupture  zone  readjust  to  equilibrium  in  a fashion  determined 
by  the  boundary  conditions  on  the  rupture  surface  once  rupture  is  initiated. 

This  readjustment  is  accomplished  by  seismic  radiation  wherever  the  medium 
is  elastic. 

Where  the  rupture  process  occurs  instantaneously,  the  elasto- 
dynamic  source  model  and  the  dislocation  source  model  are  equivalent  (Snoke 
et  al.  , 1975;  Snoke,  1975).  The  spectral  structure  obtained  using  Archambeau's 
(1968)  model  for  far-field  displacement  amplitudes  differs  from  those  obtained 
from  other  source  models,  in  that  the  spectra  obtained  from  the  elastodynamic 
model  exhibit  peaks  at  long  periods.  Snoke  (1975)  asserts  that  these  peaks  are 
actually  spurious,  and  are  due  to  physically  inadmissible  applications  of 
Archambeau's  method  of  solution. 


MOTIVATIONS  FOR  THE  DISCRIMINANTS  USED  IN  THIS  STUDY 


In  order  to  classify  seismic  events  of  unknown  origin  either 
as  earthquakes  or  as  underground  explosions,  it  is  first  necessary  to  establish 
a set  of  discriminant  measurements  on  which  the  classification  scheme  may  be 
based.  These  discriminants  must  be  chosen  to  reflect  the  differences  in  the 
source  characteristics  for  the  two  types  of  events.  Also,  because  the  dis- 
criminants are  usually  obtained  from  far-field  seismic  measurements  and 
must  be  applied  to  events  from  tectonically  different  source  regions,  the  dis- 
criminants must  be  chosen  to  be  as  insensitive  as  possible  to  variations  in- 
troduced by  differing  travel  paths  and  near-source  environments. 

Several  discriminants  have  been  developed  in  recent  years  which 
meet  the  above  requirements  to  some  degree.  Each  is  based  upon  one  (or  a 
combination)  of  various  observable  features  of  the  far-field  radiation  fields. 
Some  of  the  characteristic  differences  between  the  radiation  fields  of  the  two 
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types  of  events  which  have  been  found  to  be  useful  in  seismic  discrimination 


Differences  in  partitioning  of  radiated  energy  at  the  source 
into  shear  wave  and  compressional  wave  radiation. 


Differences  in  the  level  of  surface  wave  excitation  for  a given 


Differences  in  source  complexities  as  indicated  by  varying 
amounts  of  energy  present  in  the  P-wave  codas. 


Differences  in  the  spectral  characteristics  of  the  radiation 
fields  of  each  type  of  event. 


The  earthquake  source  models  (i.  e.  , the  dislocation  source 
model  and  Archambeau's  (1968)  elastodynamic  source  model)  mentioned  earlier 
in  this  section  may  be  described  in  terms  of  a shear  failure  along  a plane. 

In  terms  of  first  motions,  the  radiation  field  observed  from  an  earthquake 
source  is  equivalent  to  that  which  would  be  observed  by  replacing  the  earth- 
quake source  by  an  equivalent  point  double-couple  (Chinnery,  I960;  Burridge 
and  Knopoff,  1964).  This  has  also  been  verified  by  several  large-scale  ob- 
servational studies  (Hodgson  and  Stevens,  1964;  Stevens,  1969). 


On  the  other  hand,  the  underground  explosion  source  as  described 
by  Sharpe  (1942)  and  Carpenter  (1967)  may  be  represented  in  terms  of  an  im- 
pulsive, radially  symmetric  point  dilatational  source.  A significant  difference 
between  the  double-couple  source  and  the  dilatational  source  is  in  the  amount 
of  energy  partitioned  into  compressional  and  shear  wave  radiation  at  the  source. 
From  symmetry  considerations,  the  dilatational  source  corresponding  to  an 
underground  explosion  generates  only  compressional  waves  in  a homogeneous 
and  isotropic  medium  (Von  Seggern,  1972).  Shear  (SV)  bodywaves  and  Ray- 
leigh surface  waves  are  also  observed  in  the  explosion  radiation  field,  due  to 


mode  conversions  along  the  travel  path  and  the  presence  of  a free  surface 
boundary.  However,  except  for  the  effects  of  non-linear  mode  conversions 
and  near- receiver  scattering,  no  significant  horizontally  polarized  shear 
bodywave  or  Love  surface  wave  components  should  be  observed  in  the  explo- 
sion signature. 

In  contrast  to  theoretical  expectations,  significant  horizontally 
polarized  shear  wave  components  have  been  observed  in  the  radiation  fields 
of  some  explosions.  Archambeau  and  Sammis  (1970),  and  Archambeau  (1972) 
have  suggested  a mechanism  to  explain  this  phenomenon,  where  the  anomalous 
SH  component  may  be  contributed  by  a combination  of  stress  relaxation  ac- 
companying the  creation  of  a fracture  zone  in  a prestressed  medium,  and  by 
weak  zones  and  faulting  in  the  immediate  vicinity  of  the  detonation  point.  The 
latter  mechanism  leads  to  secondary  ruptures  with  definite  fault  symmetries. 
These  mechanisms  have  been  studied  by  Archambeau  and  Sammis  to  obtain 
estimates  for  the  source  parameters  of  the  Bilby  underground  explosion. 
Similar  results  have  been  obtained  by  Toksoz  et  al.  (1971),  where  the  ex- 
plosive source  is  modeled  in  terms  of  a composite  dilatational  and  double- 
couple point  source.  The  strength  of  the  double-couple  component  with  re- 
spect to  the  dilatational  contribution  is  found  to  increase  with  increasing 
strength  of  the  surrounding  medium.  Moreover,  the  orientation  of  radiation 
pattern  from  the  double-couple  has  been  found  to  be  in  general  agreement  with 
regional  faulting  trends,  which  supports  the  hypothesis  that  the  observed  SH 
wave  component  is  due  to  explosion  induced  tectonic  strain  release. 

An  automatic  spectral  fitting  procedure  which  determines  a 
minimum  least  squares  combined  dilatational  explosion  source  and  double- 
couple earthquake  source  was  developed  by  Tsai  (1972).  Since  Tsai's  method 
was  established,,  a few  studies  have  been  done  by  several  authors  (Tsai  1972a, 
1972b;  Tsai  and  Shen,  1972;  Turnbull  et  al.  , 1973,  1974,  1975;  Sun,  1976, 
1977;  Hsiao,  1978)  using  this  method.  All  of  these  studies  have  been  aimed 
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toward  the  better  understanding  and  utilization  of  the  long-period  teleseismic 
surface  waves  for  source  characterization  and  possibly  for  source  discrimi- 
nation. Turnbull  (1971,  1972)  extended  the  double-couple  source,  which  is  a 
direct  result  of  the  work  of  Ben-Menahem  and  Harkrider  (1964)  and  Harkrider 
(1970),  to  higher  order  multipolar  source.  Alexander  and  Turnbull  (1973) 
also  combined  the  fundamental  and  first  higher  mode  Rayleigh  and  Love  wave 
spectra  as  a function  of  source  parameter  variation  in  order  to  investigate  the 
possibility  of  obtaining  a more  precise  surface  wave  magnitude  measurement. 
From  the  examination  of  the  observed  Rayleigh  and  Love  wave  amplitude 
spectra  of  twenty- seven  presumed  underground  nuclear  explosions.  Sun  (1977) 
found  that  the  seismic  sources  of  those  explosions  have  various  degrees  of 
double-couple  component  in  addition  to  the  explosive  source.  Therefore,  a 
combined  source,  which  consists  of  a point  explosive  source  and  a point 
double-couple  source,  has  been  used  by  Sun  (1977)  to  model  the  explosive 
source  and  by  Hsiao  (1978)  to  model  both  explosive  source  and  earthquake 

source.  The  results  indicate  that  the  estimates  of  dip  angle  of  nine  explosions 

o o 

most  likely  occur  between  40  to  50  . This  kind  of  consistency  might  imply 
that  the  explosion  can  induce  structure  rupture  along  the  planes  of  maximum 
shear  stress  which  usually  are  oriented  around  45°.  The  results  also  show 
the  lack  of  a dilatational  component  in  long-period  explosion  surface  waves, 
random  strike  and  deeper  focus  which  suggests  that  the  double-couple  source 
of  explosions  may  be  associated  with  large  dislocations  on  a fracture  zone 
surrounding  the  source. 

Another  explanation  for  the  observed  SH  and  Love  wave  com- 
ponents in  some  explosions  is  given  by  Viecelli  (1973),  which  attributes  these 
components  to  a spalling  phenomenon.  This  is  supported  by  Von  Seggern's 
(1973)  observations  of  events  at  the  Amchitka  Test  Site.  Briefly,  Von  Seggem 
found  that  in  situ  stress  measurements  yielded  a residual  stress  level  too 
low  to  account  for  the  observed  SH  component  in  terms  of  a tectonic  release 
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mechanism.  Moreover,  delays  in  the  Rayleigh  wave  arrival  as  predicted  by 
Viecelli  (1973)  for  a spalling  source  were  observed. 


REGIONAL  PHASE  RATIO  DISCRIMINANTS 


The  theoretical  difference  between  earthquakes  and  explosions 
in  the  relative  amounts  of  energy  partitioned  into  compressional  and  shear 
wave  radiation  has  been  used  in  short-period,  regional  discrimination  studies 
by  Booker  and  Mitronovas  (1964)  and  Bell  (1978).  In  both  cases  discriminants 
were  constructed  by  taking  the  ratio  of  energy  in  a compressional  wave  veloc- 
ity window  to  the  energy  in  shear  and  short-period  surface  wave  (Lg)  velocity 
windows.  Although  the  propagation  mechanism  for  the  Lg  phase  is  not  fully 
understood,  several  explanations  have  been  advanced  which  suggest  that  the 
phase  is  related  to  shear  wave  propagation  across  a continental  structure. 


The  Lg  phase  is  a short-period  (approximately  0.  5 to  6.  0 
seconds),  large  amplitude, emergent  arrival  which  is  confined  almost  entirely 
to  continental  travel  paths,  with  an  average  source- receiver  velocity  of  about 
3.  5 km/sec  (Ewing  et  al.  , 1957;  Landers,  1978).  The  Lg  phase  is  most 
pronounced  on  the  transverse  component;  however,  it  is  also  recorded  on  the 
vertical  and  radial  components  to  some  extent.  Ewing  et  al.  (1957),  and 
Herrin  and  Richmond  (I960)  suggest  that  Lg  is  actually  a guided  shear  wave, 
traveling  in  the  upper  portion  of  the  continental  crust.  More  recent  studies 
(Knopoff  et  al.  , 1973;  Panza  and  Calcagnile,  1975;  Isacks  and  Stephens, 

1975)  suggest  that  Lg  may  be  identified  with  higher- mode  Love  wave  propaga- 
tion, which  does  not  require  a low  velocity  channel  in  the  Earth's  crust.  Lan- 
ders (1978)  attributes  the  vertical  part  of  Lg  to  higher-mode  Rayleigh  waves 
and  scattering  near  the  receiver.  In  each  case,  however,  the  excitation  of 
Lg  is  directly  coupled  to  the  partitioning  of  energy  into  shear  wave  radiation 
at  the  source,  suggesting  its  usefulness  in  a seismic  discrimination  effort. 


Of  the  various  phase  energy  ratios  studied  by  Booker  and 
Mitronovas  (1964)  and  Bell  (1978),  the  most  powerful  (in  a seismic  discrimina 
tion  context)  were  found  to  be: 


The  ratio  of  energy  in  the  Pg  window  (v  (average  velocity) 
6.  9 km/sec  - 4.  9 km/sec)  to  the  energy  in  the  Lg  window 
(v  = 3.  6 km/ sec  - 3.  2 km/ sec). 


The  ratio  of  energy  in  the  Pg  window  to  the  energy  in  a shear 
bodywave  window  (4.9  km/ sec  -3.6  km/ sec). 


The  ratio  of  energy  in  a window  containing  Pn  and  Pg  (v  = 
first  arrival  - 4.  6 km/sec)  to  the  energy  in  a window  contain 
ing  shear  bodywaves  and  Lg  (v  = 4.6  km/sec  -2.0  km/sec) 


These  discriminants  were  applied  to  a set  of  20  North  American 
earthquakes  and  27  Nevada  Test  Site  underground  explosions.  Threshold  F- 
statistics  were  obtained  which  support  the  hypothesis  that  the  earthquake  and 
explosion  group  means  for  these  discriminants  were  unequal  at  the  99  per- 
cent confidence  level  in  all  three  cases  (Bell,  1978). 


Where  short-period  S and  Lg  data  are  available,  regional  dis 
criminants  based  upon  the  difference  in  energy  partitioning  are  used.  How- 
ever, instead  of  using  ratios  of  total  energy  as  did  previous  studies,  the  dis 
criminants  are  constructed  by  taking  the  difference  between  magnitude  mea- 
surements for  the  P,  S,  and  Lg  phases.  The  phase  magnitudes  are  based 
upon  broadband  amplitude-period  ratio  measurements,  and  a:  <*  scaled  to  m 

b 

according  to  the  following  magnitude  formulas: 


3.  82  + 2.  00  logA  + log(A/T)  (Evernden,  1967)  (II- 1) 


(Fitch  et  al.  , 1978) 


mgn  = 3.  79  + 1.  80  log  A + log  (A/T) 

mT  = 3.  30  + 1.66  logA  + log  (A/T)  (Nuttli,  1973) 

J-’g 


After  evaluating  the  magnitudes,  the  discriminants  Am  and  Am  are 

* c* 

formed,  where 


Am  = rn  - m _ 

1 Sn  Pn 

Am  = rn  _ rn 

2 Lg  Pn 


It  may  be  seen  that  these  discriminants  are  analogous,  in  a formal  sense,  to 

the  teleseismic  M - m,  discriminant, 
s b 

The  magnitude  formulas  for  Pn,  Sn  and  Lg  are  based  upon 
the  observed  attenuation  characteristics  of  Eastern  North  America,  and  have 
also  been  found  to  be  applicable  to  many  Eurasian  events  (Nuttli,  1973;  Fitch 
et  al.  , 1978).  Because  of  the  complex  structure  of  tectonically  active  regions 
such  as  those  found  in  the  Western  United  States  and  Kuriles -Kamchatka,  the 
magnitude  formulas  may  be  unreliable  when  applied  to  events  from  these  areas. 
Magnitude  formulas  which  apply  to  tectonically  active  and  rift  zones  could  be 
developed  for  a properly  regionalized  data  base.  However,  the  event  data  base 
used  in  the  present  study  contains  too  few  events  to  adequately  represent  each 
source  region,  thus  it  was  felt  that  greater  harm  would  be  done  by  attempting 
to  regionalize  this  data  base  than  by  considering  it  in  toto. 

D.  THE  TELESEISMIC  Mg-  rr^  DISCRIMINANT 

In  addition  to  being  more  efficient  in  generating  shear  body- 

waves,  the  shearing  source  corresponding  to  an  earthquake  is  also  a more 

efficient  generator  of  surface  waves  than  the  dilatational  explosive  source 

(Douglas  et  al.  , 1971;  Gilbert,  1973).  This  forms  the  basis  for  the  M - m 

s b 
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discriminant  introduced  by  Press  et  &1.  (1963),  which  has  proven  to  be  one  of 
the  most  effective  teleseismic  discriminants  available. 

In  addition  to  measuring  the  difference  in  the  relative  excita- 
tion of  bodywave  and  surface  wave  radiation,  the  discriminant  also 

reflects  differences  in  the  relative  amounts  of  high  and  low  frequency  energy 
generated  by  the  seismic  source.  A general  feature  of  theoretical  far-field 
displacement  spectra  for  explosive,  shear  dislocation  and  stress  relaxation 
sources  is  the  characteristic  roll-off  in  the  spectra  of  bodywaves  and  surface 
waves  at  frequencies  above  some  critical  corner  frequency  which  varies  in- 
versely with  the  source  dimension  (Hanks  and  Wyss,  1972;  Aki,  1967). 

Above  the  corner  frequency,  the  displacement  spectrum  must 
decrease  at  least  as  quickly  as  co  ^ to  ensure  bounding  of  the  energy  inte- 
gral. Aki  (1967)  discusses  two  statistical  source  models,  one  where  the  am- 
plitude spectrum  decays  as  w ^ at  high  frequencies,  and  another  due  to 

_ 3 

Haskell  (1966),  which  decays  as  w . Geller  (1976)  has  also  proposed  a 

model  where  fault  width,  as  well  as  length,  and  the  direction  of  propagation 

are  explicitly  taken  into  account,  resulting  in  different  average  spectra  for 

bodywaves  and  surface  waves.  Above  the  comer  frequency  for  width  the 

_3 

average  spectrum  decays  as  u>  .In  the  region  between  the  corner  fre- 

-2 

quencies  for  length  and  width,  the  spectrum  begins  to  decay  as  w , and  be- 
low the  corner  frequency  for  length,  the  average  spectrum  is  essentially  flat. 


For  each  of  the  representations  described  above,  the  corner 

frequencies  decrease  with  increasing  magnitude.  This  tends  to  segment  the 

theoretical  Mg  versus  m^  line,  where  the  slope  of  the  line  varies  according 

to  whether  M , m,  or  both  M and  m,  are  measured  above  or  below  the  the- 
s b s b 2 

oretical  corner  frequency.  For  > 6.  5,  the  slope  found  by  Aki  for  o» 
model  closely  follows  the  Gutenberg  and  Richter  (1956)  Mg  versus  m^  line, 
but  diverges  for  Mg  < 6.  5. 


II- 11 


I 


r 


t ! 

K 


Liebermann  and  Pomeroy  (1969)  utilized  Aki's  (1967)  theoreti- 
cal scaling  relations  for  earthquakes,  and  a scaling  relation  for  explosions 
where  magnitude  scales  as  the  common  logarithm  of  the  cavity  radius.  For 
Mg  and  measurements  at  20  seconds  and  1 second  period,  respectively, 
the  dependence  of  Mg  and  upon  source  dimension  for  earthquakes  is  given 
by  Liebermann  and  Pomeroy  (1969): 


/ log  L + constant 

Ms  > 6.  5 

Mg  — < 2 log  L + constant 

3.  5 < Mg  < 6.  5 

\ 3 log  L + constant 

Mg  < 3.5 

£ ( log  L + constant 

Mg  > 3.  0 

I mb  ~ i 

( 2 log  L + constant 

Ms<  3.  0 

where 

L is  the  fault  length. 
Similarly,  for  explosions, 

Ms  ~ 3 log  a + constant 


m,  ~ 3 log  a + constant 

D 


% 


< 6.  0 


where 


a is  the  radius  of  the  effective  cavity. 


The  results  for  M and  m,  may  be  combined  to  give: 
s b 

Earthquakes: 


Mg  = m^  + constant 

M > 

s 

o 

in 

Mg  = 2m_b  + constant 

3.  5 < M < 
s 

6.  5 

M = 1.5m,  + constant 

M < 

3.  5 

s b 

s 

Explosions: 

M = m,  + constant 
s b 

m < 
b 

6.  0 . 
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Thus,  for  M > 3.  5 the  earthquake  and  explosion  populations  diverge,  and  for 
s 

M < 3.  5,  the  two  populations  appear  to  merge, 
s 

These  theoretical  results  indicate  that  larger  earthquakes  tend 
to  excite  lower  frequencies  more  efficiently  than  higher  frequencies,  thus  en- 
hancing the  M^-  discriminant  for  3.  5 < Mg<  6.  5.  Clearly  .additional  data 
for  low  magnitude  events  are  required  to  establish  the  behavior  of  m^  at 
low  magnitudes  and  the  validity  of  Aki's  theoretical  scaling  law. 

Many  studies  have  been  conducted  to  evaluate  the  effectiveness 
of  the  Mg-  discriminant  involving  numerous  events  from  many  different 
source  regions  (Thirlaway,  1968;  SIPRI,  1968;  Basham,  1969;  Liebermann 
and  Pomeroy,  1969;  Capon  et  al.  , 1969;  Evernden,  1969;  Evernden  and  Fil- 
son,  1971;  Marshall  and  Basham,  1972;  Peppin  and  McEvilly,  1974;  and 
Evernden,  1975;  to  name  just  a few).  M^-  m^  has  been  found  to  reliably  sep- 
arate earthquakes  and  underground  explosions  for  m^  values  as  low  as  4.  75 
(Evernden,  1969)  to  5.00  (Thirlaway,  1968;  SIPRI.  1968). 

Estimates  of  bodywave  and  surface  wave  magnitudes  may  be 
influenced  by  factors  which  are  not  specifically  source  related.  These  factors 
tend  to  introduce  scatter  and  bias  into  the  magnitude  measurements  and  tend 
to  reduce  the  effectiveness  of  Mg-  in  the  context  of  a world-wide  discrim- 
ination scheme.  Considerable  research  has  been  directed  toward  identifying 
sources  of  magnitude  bias  and  removing  or  minimizing  this  bias  where  pos- 
sible. 

One  source  of  error  occurs  where  differing  attenuation  char- 
acteristics for  various  source- receiver  paths  are  not  accounted  for.  Con- 
sistent magnitude  differences  have  been  found,  for  example,  for  teleseismic 
events  measured  in  the  western  United  States  and  for  those  measured  in  the 
eastern  United  States  (Evernden  and  Clark,  1970;  Booth  et  al.  , 1974;  Der, 
1977).  Liebermann  and  Pomeroy  (1969)  and  Basham  (1969)  also  note 
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systematic  differences  in  the  Mg-  nr^  relations  (i.  e.  , the  least  squares 
Mg-  lines)  for  events  in  the  western  United  States  as  compared  to  those  in 
Aleution- Kamchatka  area,  southern  Algeria,  central  Asia,  and  Novaya  Zemlya. 
For  these  regions  the  observed  shift  in  the  Mg-  nr^  lines  has  been  found  to 
range  from  0.  3 to  0.  5 magnitude  units  (Evernden  and  Filson,  1971;  Solomon, 
1972). 

Regional  dependence  of  Mg-  m^  has  been  explained  in  terms  of 
differing  levels  of  anelastic  attenuation  associated  with  each  of  the  source  re- 
gions under  consideration.  This  suggests  the  need  for  consideration  of 

on  a regional,  rather  than  global,  basis  (Liebermann  and  Pomeroy, 
1969;  Evernden  and  Filson,  1971;  Marshall  and  Basham,  1972).  In  addition, 
Der  (1977)  suggests  that  a world-wide  map  of  Q (or  t*)  values  be  compiled 
which  would  allow  for  correction  of  Q variations  under  both  source  and  re- 
ceiver. 

Source  depth  also  influences  surface  wave  excitation,  as  noted 
by  Evernden  (1975).  In  general,  deeper  earthquakes  result  in  less  surface 
wave  excitation  relative  to  those  occuring  at  shallow  depths.  Several  methods 
have  been  studied  for  estimating  source  depth,  such  as  the  application  of  cep- 
stral  analysis  and  matched  filter  techniques  to  depth  phases  (Page,  1976)  and 
by  fitting  source  models  to  surface  waves  (Sun,  1976,  1977;  ahd  Hsiao,  1978). 
Such  techniques  might  be  employed  in  a screening  process  (Basham  and  Anglin, 
1973;  Bell,  1978)  to  quickly  identify  events  as  natural  earthquakes  which  occur 
at  depths  greater  than  a few  tens  of  kilometers,  or  which  occur  at  sea  or  other 
areas  where  bomb  emplacement  is  highly  impractical. 

E.  COMPLEXITY 

THe  underlying  concept  of  the  teleseismic  P wave  complexity 
discriminant  is  that  the  explosive  source  may  be  represented  as  a transient 
disturbance  with  a very  short  rise  time  (Liebermann  and  Pomeroy,  1969). 
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As  a result,  the  first  arrival  consists  of  a short-period  P wave  with  little 
trailing  coda.  In  comparison,  an  earthquake  is  also  represented  as  a trans- 
ient phenomenon;  however,  the  earthquake  has  a finite  time  duration  as  the 
disturbance  propagates  along  the  fault.  This  finite  time  duration  tends  to 
produce  a significant  coda  following  the  initial  P wave  arrival.  An  example 
of  this  is  shown  in  Figure  II- 1,  where  the  upper  trace  is  an  eastern  Kazakh 
underground  nuclear  explosion  and  the  lower  trace  is  an  earthquake  from  the 
same  region. 

The  complexity  discriminant  consists  of  measurements  of  the 
energy  in  the  P wave  coda  taken  at  teleseismic  distances.  Complexity  mea- 
surements are  generally  not  made  at  regional  distances  because  the  presence 
of  closely  arriving  crustal  phases  (i.  e.  , Pn,  Pg,  and  P*)  make  the  signatures 
of  both  earthquakes  and  explosions  appear  complex.  The  conventional  com- 
plexity measurement  is  the  inverse  ratio  of  the  P wave  energy  in  the  first  5 
seconds  after  the  initial  P wave  arrival  to  the  energy  in  the  following  30  sec- 
onds (Evernden,  1969;  Lambert  et  al.  , 1970).  These  time  windows  are  not 
invariant,  as  demonstrated  by  Anglin  (1971),  who  found  that  gates  from  0 to 
2 seconds  and  2 to  35  seconds  were  optimum  for  measuring  complexity  on 
Eurasian  data  as  recorded  at  the  Yellowknife  array  (YKA)  in  Canada. 

F.  VARIABLE  FREQUENCY  MAGNITUDE  METHODS 

Variable  frequency  magnitude  methods  have  been  applied  to 
short-period  teleseismic  P waves  to  discriminate  explosions  from  earthquakes. 
Bache  et  al.  (1975)  thoroughly  evaluated  a short-period  variable  frequency 
magnitude  discriminant  described  by  Savino  and  Archambeau  (1974),  Bache 
et  al.  (1974),  and  Archambeau  et  al.  (1974).  Their  results  indicated  operat- 
ing characteristics  similar  to  the  M versus  m discriminant.  Their  method 


FOR  EASTERN  KAZAKH  EXPLOSION  AND  EARTHQUAKE 
AS  RECORDED  AT  MAIO 


trend  line  is  drawn  through  the  earthquake  population,  then  the  distance  of 
explosion  coordinates  from  the  earthquake  trend  line  is  approximately  three 
quarters  of  a magnitude  unit.  The  explosions  are  well  separated  from  earth- 
quakes by  enhanced  high  frequency  energy  and  by  less  low  frequency  energy. 
The  effect  of  both  of  these  changes  is  more  normal  displacement  of  the  ex- 
plosion points  from  the  earthquake  trendline. 

Sax  (1976)  extended  the  variable  frequency  magnitude  method 
to  multivariate  analysis  of  low  and  high  frequency  magnitude  measurements 
and  to  the  measured  slope  of  the  magnitude  spectrum.  The  measurements  of 
magnitudes  at  each  frequency  were  reduced  to  a single  variable  which  is  in- 
dependent of  event  size.  This  was  done  by  subtracting  from  each  high  or  low 
frequency  magnitude  measurement  the  average  magnitude  in  the  band  between 
1 Hz  and  1. 75  Hz.  The  resultant  low  frequency  and  high  frequency  magnitude 
and  roll-off  measurements  were  treated  as  single  ten-dimensional  vector 
measurements.  A scalar  mean  and  standard  deviation  was  computed  for  each 
element  of  the  vector  to  represent  a 'normal'  population  of  shallow  earthquakes 
To  compute  a discriminant  'z'  vector,  each  unknown  event  vector  measure- 
ment of  magnitude  and  roll-off  was  reduced  to  a unit  normal  statistic  by  sub- 
tracting the  mean  and  dividing  by  the  standard  deviation  of  the  'normal'  earth- 
quake population.  The  final  step  was  to  determine  one  or  more  mean  z- 
discriminant  vectors  which  point  to  explosion  populations.  In  that  way,  each 
unknown  z- discriminant  vector  is  projected  in  the  direction  of  vectors  point- 
ing toward  known  explosion  populations.  The  results  indicated  excellent 
separation  of  NTS  presumed  explosions,  eastern  Kazakh  presumed  explosions, 
and  shallow  earthquakes.  The  vectors  pointing  toward  NTS  and  eastern  Ka- 
zakh presumed  explosions  were  nearly  orthogonal. 


c. 


INSTANTANEOUS  PHASE  AND  FREQUENCY  MEASUREMENTS 


The  theoretical  result  that  earthquakes  tend,  in  general,  to 
excite  lower  frequency  radiation  more  efficiently  than  explosions  from  the 
same  source  region  forms  the  basis  for  the  mean  instantaneous  frequency 
discriminant.  The  mean  instantaneous  frequency  has  been  studied  in  conjunc- 
tion with  the  mean  phase  standard  deviation  by  Unger  (1978a).  Used  as  a two- 
variate  discriminant  these  parameters  were  found  to  significantly  separate 
NORSAR  single  site  recordings  of  Eurasian  earthquakes  and  Russian  presumed 
underground  explosions. 

It  has  been  demonstrated  that  the  presence  of  early  secondary 
arrivals  such  as  pP  can  be  determined  from  rapid  changes  in  the  time  series 
of  instantaneous  amplitude,  phase,  and  frequency  (Farnbach,  1975;  Unger, 
1976).  The  number  of  secondary  arrivals  may  be  quantized  by  measuring  a 
regression  analysis  residue  of  the  instantaneous  phase  time  series,  as  de- 
scribed in  Section  HI  of  this  report,  and  by  Unger  (1978a).  This  residue  is 
termed  the  phase  standard  deviation. 

By  plotting  a measure  of  mean  instantaneous  frequency  against 
the  mean  phase  standard  deviation  for  the  first  two  seconds  immediately  fol- 
lowing the  P wave  onset,  Unger  (1978a)  found  significant  separation  between 
populations  of  Russian  presumed  underground  explosions,  Eurasian  earth- 
quakes, and  Nevada  Test  Site  (NTS)  underground  explosions.  Relative  to  the 
Eurasian  earthquakes,  the  Russian  presumed  explosions  exhibited  both  a high- 
er dominant  frequency  and  a greater  number  of  early  secondary  signals  as 
indicated  by  a larger  mean  phase  standard  deviation.  In  contrast,  the  NTS 
events  showed  near-monochromatic,  low  frequency  waveforms,  suggesting 
a lack  of  clearly  defined  early  multiple  signals. 

The  mechanism  whereby  Russian  presumed  explosions  might 
exhibit  more  early  multiple  arrivals  than  Eurasian  earthquakes  or  NTS  ex- 
plosions is  not  well  understood.  However,  for  the  same  data  base,  identical 
event  classifications  were  independently  obtained  by  Sax  (1976). 
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SECTION  m 

MEASUREMENT  OF  EVENT  IDENTIFICATION  PARAMETERS 

A.  OVERVIEW  OF  THE  MEASUREMENT  PROGRAM 

The  purpose  of  this  section  is  to  describe  the  manner  in  which 
the  various  event  identification  parameters  are  measured  on  the  data  record- 
ed at  the  individual  stations.  The  program  which  performs  these  measure- 
ments, called  TISSPROG,  is  described  by  Schmidt  and  Wilson  (1978).  This 
discussion  will  concentrate  on  those  aspects  of  the  program  which  deal  with 
the  event  identification  problem. 

The  concept  underlying  the  creation  of  this  program  is  that  it 
is  more  efficient  to  make  one  computer  run  which  will  measure  all  desirea 
event  identification  parameters  than  to  run  a series  of  programs  each  of 
which  measures  one  parameter  or  a group  of  related  parameters.  The  one- 
run  approach  minimizes  the  amount  of  time  required  to  prepare  for  the  com- 
puter run  as  well  as  the  actual  amount  of  computer  time  used.  As  an  addi- 
tional benefit,  the  one- run  approach  produces  an  output  containing  all  desired 
measurements.  Thus,  the  analyst  is  not  required  to  spend  time  assembling 
the  results  of  multiple  computer  runs  before  proceeding  with  his  analysis. 

In  the  most  general  sense,  the  operation  of  the  program  is  to  compute  neces- 
sary event  parameters  such  as  travel  times,  determine  the  type  of  measure- 
ments to  be  made  for  each  input  waveform,  make  these  measurements,  save 
the  results  in  a manner  suitable  to  the  needs  of  the  analyst,  and  repeat  the 
procedure  for  the  next  selected  event. 

In  order  to  minimize  the  amount  of  information  to  be  input  by 
the  analyst, a number  of  parameters  have  been  built  into  the  program.  These 
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parameters  are: 

• Station  coordinates 

• Station  names  (four-letter  designators) 

• Teleseismic  P,  S,  LQ,  and  LR  travel  time  equations 

• Station  tectonic  classification  - active,  inactive,  or  rift 

• Short-period  noise  magnitudes  for  the  SRO  and  ASRO  stations 

• Short-period  and  long-period  detection  thresholds 

• Veith  and  Clawson  (1972)  P factors  for  bodywave  magnitude 

computation 

• System  response  corrections 

• Corner  frequencies  of  short-period  and  long-period  filters 

• Empirical  dispersion  curve. 

The  purpose  of  each  of  these  built-in  parameters  will  be  de- 
scribed in  the  course  of  this  section. 

The  generalized  flow  of  the  event  identification  parameter 
measurement  program  is  shown  in  Figure  111- 1 . The  basic  parameters  input 
to  the  program  are  date,  origin  time,  and  coordinates  of  the  event  (input 
parameter  card)  and  station  number,  seismogram  number,  data  type  (SP  or 
LP),  channel  number,  type  of  waveform  (P,  S,  .or  LR),  and  resample  rate 
(seismogram  cards),  where  one  seismogram  card  is  required  for  each  wave- 
form edit  on  the  Area  of  Interest  (AI)  digital  magnetic  tape. 

After  opening  and  positioning  the  input  AI  tape  and  the  output 
event  waveform  tape  and  discriminant  tape,  the  program  reads  the  event 
parameter  card  and  the  first  seismogram  card.  (Other  parameters  such  as 
depth  and  bodywave  magnitude,  if  available,  can  be  entered  via  the  event 
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parameter  card.  However,  the  only  event  parameters  necessary  for  the  event 
identification  study  are  date,  origin  time,  and  event  coordinates.  ) The  pro- 
gram now  determines  the  station  coordinates  from  the  station  number  and 
built-in  station  coordinates.  Using  this  information  and  the  read-in  event 
coordinates,  station-to-event  azimuth  and  distance  (delta)  are  next  calculated. 
The  built-in  travel  time  equations,  input  origin  time,  and  computed  delta  are 
then  used  to  compute  the  arrival  times  of  the  propagation  modes  at  the  station. 
A Rayleigh  wave  end  time  based  on  a velocity  of  2.  5 km/sec  and  the  Rayleigh 
wave  arrival  time  are  also  computed.  With  this  information,  a long-period 
edit  start  time  is  computed  from 
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where  J is  a decimation  index  (original  sample  rate/desired  sample  rate). 
The  maximum  length  of  the  edit  is  set  at  2048  points.  Similar  edit  times  are 
computed  for  the  short-period  phases.  With  this  information,  the  program 
reads  the  data  from  AI  tapes  as  specified  by  the  seismogram  number.  Checks 
are  made  and  edit  start  and  end  points  are  reset  as  required  by  the  actual 
amount  of  data  recorded  on  the  AI  tape  for  that  seismogram  number.  In  the 
case  of  array  data,  data  recorded  at  only  one  site  is  read  into  the  program. 

At  all  arrays  except  ILPA,  the  site  used  is  the  site  nearest  the  center  of  the 
array.  At  ILPA,  site  7 was  used,  since  this  is  the  only  site  of  the  array  re- 
cording both  short-period  and  long-period  data. 

The  next  step  is  to  compute  and  remove  the  mean  from  each 
data  trace  and,  in  the  case  of  three  component  long-period  data, to  rotate  the 
data  from  its  recorded  vertical,  north,  east  configuration  to  a vertical,  trans- 
verse, radial  configuration  as  determined  by  the  computed  azimuth.  At  this 
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point,  the  data  are  ready  for  the  measurement  stages  of  the  program.  The 
program  therefore  tests  the  input  data  type  flag  and  directs  the  program  flow 
toward  the  short-period  or  long-period  data  measuring  algorithms  as  appro- 
priate. 

B.  SHORT- PERIOD  EVENT  MEASUREMENT 

The  short-period  discriminant  measurement  phase  of  TISSPROG 
proceeds  as  illustrated  in  Figure  EII-2.  The  short-period  discriminant  mea- 
surements are  directed  by  a supervisory  routine  (subroutine  SPEED),  which 
directs  control  to  each  of  the  various  discriminant  measurement  and  data 
manipulation  routines  as  required.  In  the  case  of  magnitude  measurements, 
only  log  (A/T)  or  log  A values  are  calculated  in  SPEED.  These  values  are 
then  placed  in  the  event  header  for  magnitude  calculation  in  subroutine  FINISH. 

1.  Signal  Processing  Prior  to  Discriminant  Measurements:  Filter- 

ing, Signal  Detection,  Instrument  Response  Removal,  and  Edit 
Compression 

The  first  step  in  the  short-period  discriminant  measurement 
phase  is  to  bandpass  filter  the  input  raw  data  trace  over  a broadband  frequency 
range.  This  tends  to  subdue  high  frequency  spikes  (i.  e.  , artifacts  due  to  the 
instrumentation)  and  high  frequency  components  in  the  data  to  avoid  aliasing. 

Because  the  frequency  interval  is  fairly  broadband,  a rectangular  boxcar  fil- 
ter is  used  without  risk  of  severe  ringing  problems.  In  the  case  of  data  sam- 
pled at  20  Hz,  the  filter  cutoffs  are  0.  3 Hz  to  9.  9 Hz; where  the  data  are  sam- 
pled at  10  Hz,  the  filter  cutoffs  are  set  at  0.  3 Hz  and  4.  9 Hz. 

After  the  raw  data  trace  has  been  filtered,  the  next  step  con- 
sists of  running  the  automatic  short-period  signal  detector  described  by 
Unger  (1978a).  The  procedure  for  detecting  and  timing  the  onset  of  short- 
period  signals  is  illustrated  in  Figure  HI- 3.  Over  a previously  specified 
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warm-up  time  (10  seconds  for  events  in  this  study),  the  peak  noise  envelope 
is  established.  The  peak  noise  envelope  is  then  cosine  tapered  over  subse- 
quent points.  An  envelope  value  exceeding  the  tapered  peak  results  in  an  up- 
date in  the  noise  peak.  Where  a signal  is  declared,  no  noise  peak  update  take 
place  until  the  signal  has  been  declared  to  be  terminated. 

A declared  signal  detection  results  whenever  a trigger  proba- 
bility exceeds  a specified  threshold  (TH1  in  Figure  III- 3,  set  to  0.  3 in  this 
study).  The  trigger  probability  is  the  probability  that  the  envelope  in  a four- 
second  leading  time  window  is  greater  than  the  tapered  peak  noise  envelope. 
When  this  probability  reaches  a maximum  the  algorithm  looks  for  the  first 
signal  envelope  peak.  Where  the  ratio  of  the  first  signal  envelope  peak  to  the 
tapered  noise  envelope  peak  exceeds  the  signal-to-noise  threshold  TH2  (set  to 
6.  0 dB),  the  signal  detection  is  confirmed  and  a frequency  dependent  stepback 
procedure  is  performed  to  determine  the  signal  onset  time. 

The  stepback  procedure  used  in  timing  the  signal  onset  is 
shown  in  Figure  HI-4.  It  is  based  upon  the  observation  that  in  most  cases  the 
first  signal  envelope  peak  (t^)  occurs  within  one  signal  period,  typically  at 
about  3/4  period  after  the  signal  onset  (tQ).  For  waveforms  of  high  signal- 
to-noise  ratio,  the  onset  time  is  found  by  detecting  the  first  maximum  or 
minimum  of  the  signal  waveform  at  time  (t3).  This  is  followed  by  a stepback 
of  1/4  period.  For  low  signal-to-noise  waveforms,  however,  the  first  quar- 
ter period  of  the  signal  may  be  obscured  by  noise.  Where  this  is  the  case, 
the  algorithm  steps  back  3/4  mean  period  from- the  first  signal  envelope  peak 
(t4)  to  the  point  1 2 at  0.  8 mean  period  before  The  search  for  the  first 
quarter  period  peak  or  trough  of  the  seismic  signal  then  begins  at  t_  and  ends 
as  a detected  first  quarter  period  when  its  maximum  or  minimum  exceeds  a 
threshold.  In  this  study,  this  third  threshold,  TH3,  is  set  at  one  dB  above 
the  peak  noise  in  a one  second  gate  (tj.t^).  In  the  case  where  the  signal-to- 
noise  threshold  (TH2)  is  not  satisfied,  the  detection  is  abandoned,  and  the 
noise  peak  value  is  updated  with  the  current  peak. 
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The  signal  end  time  is  taken  as  the  moment  of  the  first  envelope 
minimum  occurring  either  after  the  probability  of  the  signal  peak  exceeding 
the  noise  peak  falls  below  threshold  TH1,  or  after  the  signal  duration  exceeds 
a specified  maximum,  whichever  is  first.  If  this  envelope  minimum  is  great- 
er than  the  tapered  noise  peak,  the  noise  peak  envelope  is  updated,  and  the 
signal  detection  process  resumes  as  before. 

In  the  evaluation  by  Unger  (1978a)  it  was  determined  that  the 
optimal  warm-up  time  for  this  detector  is  on  the  order  of  20  to  30  seconds, 
and  that  the  signal-to-noise  detection  threshold  is  optimal  (in  terms  of  the 
detection  ratio-false-alarm  rate  tradeoffs)  at  about  3 dB.  However,  the  con- 
struction of  the  data  base  used  in  this  study  requires  that  the  detector  warm- 
up time  be  uniformly  reduced  to  10  seconds,  as  phase  arrivals  placed  at  the 
beginning  of  the  edits  are  not  unusual.  A possible  inflation  of  the  false-alarm 
rate  accompanying  the  reduction  of  the  warm-up  time  is  compensate  ! by  a 
corresponding  increase  of  the  signal-to-noise  threshold  from  3 dB  to  6 dB. 

As  a result,  the  detector  operates  suboptimally;  however,  this  may  be  cor- 
rected by  a more  judicious  choice  of  starting  times  in  the  original  event  edit- 
ing process. 

After  the  signal  detector  has  been  run,  the  instrument  response 
is  removed  as  follows.  First,  the  corner  frequencies  of  system  response 
poles  are  determined  from  the  system  response  curve  by  graphically  finding 
the  frequency  of  intersection  of  straight  lines  representing  the  18,  12,  6,  0, 

-6,  -12,  and  -18  dB/octave  slopes  of  the  system  response  curve.  Each  of 
these  poles  are  removed  by  applying  a first  order  difference  de-pole  operator 
to  the  data  to  remove  each  pole.  The  result  is  scaled  to  a normalized  re- 
sponse of  one  at  the  calibration  frequency  (1  Hz  for  short-period,  0.04  Hz  for 
long-period).  The  effect  of  the  de-pole  operation  is  to  straighten  out  the  re- 
sponse curve  such  that  the  system  response  is  flat  with  respect  to  ground 
motion.  An  18  dB/octave  slope  of  the  ground  displacement  corresponds  to  a 


flat  response  to  the  first  derivative  of  acceleration;  a 12  dB /octave  slope,  to 
acceleration;  a 6 dB /octave  slope,  to  velocity;  and  a flat  or  0 dB /octave  slope, 
to  displacement.  Next,  a set  of  narrowband  filters  are  used  to  integrate  the 
de-poled  data  in  such  a way  as  to  yield  a flat  ground  displacement  response. 
Prior  to  application  of  the  narrowband  filters  the  de-poled  data  are  bandpass 
filtered.  For  short-period  data,  the  filter  limits  are  between  0.  16  Hz,  the 
low  -6  dB  point  of  the  first  narrowband  filter  (centered  at  0.  3 Hz)  and  a fre- 
quency just  below  Nyquist  frequency.  For  long-period  data  the  limits  are 
between  0.  010  and  0. 117  Hz.  Integration  of  the  data  trace  is  approximated 
by  narrowband  filtering  the  data  and  summing  the  scaled  results.  The  num- 
ber of  filters  used  depends  on  which  poles  have  been  removed  from  the  data, 
that  is  on  the  slope  of  the  de-poled  data.  For  example,  if  the  de-poled  data 
has  a slope  of  18  dB/octave,  which  corresponds  to  the  first  derivative  of  ac- 
celeration, three  integrations,  therefore,  three  filters,  would  be  necessary 
in  order  to  correct  the  data  to  ground  displacement.  For  short-period  data. 


the  first  and  lowest  filter  is  centered  at  a log^  frequency  of  -0.  7 and  subse- 


quent filters  are  spaced  at  increasing  logjQ  frequencies  in  0.  2 increments. 


The  effective  bandwidth,  the  width  between  the  -6  dB  points,  of  the  filters  are 
0.  452  times  the  center  frequency  which  allows  the  -6  dB  points  of  the  filters 
to  just  overlap.  The  filtered  data  are  scaled  to  be  one  at  1 Hz.  For  long- 
period  data  the  filters  start  with  a logj^  center  frequency  of  -1. 70  and  subse- 
quent filters  are  centered  at  even  intervals  of  0.  10  in  log  frequency.  The  ef- 
fective bandwidth  is  0.  222  times  the  center  frequency.  This  again  allows  the 
-6  dB  points  to  overlap.  The  filtered  data  are  scaled  to  be  one  at  0.  04  Hz. 

In  both  cases,  the  resultant  scaled  output  for  all  the  filters  is  then  summed 
and  added  to  the  broadband  filtered  data  which  has  had  poles  put  in  at  each 
of  the  high  -6  dB  points  of  the  filters.  A pole,  1/  | l-z|  , is  approximated 
by  1/|  l-az|  where  z is  the  delay  operator  e lW*  and  a is  the  reflection 
coefficient.  The  poling  operation  is  also  scaled  to  be  one  at  the  previously 
stated  frequencies  for  short-  and  long-period  data. 


In  summary,  the  first  order  real  poles  and  zeros  used  to  rep- 
resent the  system  response  by  fitting  the  amplitude  response  curve  are  re- 
moved by  stable  inverse  operators  with  a normalized  response  of  one  at  the 

calibration  frequencies.  To  avoid  round-off  errors  in  the  inverse  of  the 
N 

zeros  (1-z)  (6N  dB/octave  response),  the  low  frequency  response  is  synthe- 

sized by  adding  appropriately  scaled  narrowband  filters  which  were  demon- 
strated to  pass  the  signal  without  seriously  distorting  the  amplitude  of  signals. 

Frequencies  above  the  bank  of  low  frequency  filters  are  flattened  by  inverting 

N . -N 

(1  - a z)  which  has  a stable  inverse.  The  corner  of  the  inverse  (1-arz)  was 

designed  to  coincide  with  the  corner  of  the  highest  frequency  narrowband  fil- 
ter. The  total  inverse  was  obtained  by  adding  the  high-pass  inverse  output 
to  the  properly  scaled  narrowband  filter  outputs.  The  total  inverse  was  de- 
signed to  be  accurate  to  approximately  0.  1 magnitude  units.  This  method 
was  designed  to  avoid  the  large  drift-errors  commonly  associated  with  the 
N-fold  (N  usually  3)  integrations  of  the  data  required  to  derive  ground  uis 
placements.  Although  the  method  appeared  to  work  satisfactorily  with  the  few 
examples  examined,  more  data  needs  to  be  examined  to  verify  this  method 
since  it  is  a critical  step  in  obtaining  accurate  frequency  dependent  magnitude 
measurements. 

In  order  to  avoid  storage  of  unnecessary  amounts  of  noise  data 
on  the  output  event  tapes  and  to  facilitate  analyst  checking  of  automatic  de- 
tections, it  was  found  desirable  to  compress  the  short-period  edit  lengths. 

In  the  case  of  P and  S wave  edits,  the  edit  length  is  shortened  to  51. 2 seconds. 
For  those  events  where  the  short-period  broadband  detector  is  triggered,  the 
edit  is  centered  on  the  detected  phase  arrival  time.  For  non-detected  events, 
the  edit  is  centered  about  the  expected  phase  arrival  time,  determined  by  a 
linear  fit  to  the  Jeffreys- Bullen  travel  time  curves  in  the  distance  range  of 
interest.  In  the  case  of  Lg  edits,  where  the  short-period  broadband  detector 
does  not  trigger,  the  edit  is  centered  about  the  point  corresponding  to  an 


average  source- receiver  velocity  of  3.  6 km/sec.  Edit  compressions  for  both 
detections  and  non- detections  are  illustrated  schematically  in  Figure  III- 5- 

During  the  edit  compression,  several  checks  are  performed 
to  ensure  that  the  proposed  compressed  edit  is  at  least  partially  contained 
within  the  original  edit.  This  allows  for  a limited  amount  of  recovery  in  the 
event  that  insufficient  data  are  supplied  on  the  original  edit  tapes  by  further 
shortening  the  compressed  edit.  For  events  where  the  broadband  detector  is 
triggered,  at  least  10  seconds  of  data  must  be  available  before,  and  25.  6 sec- 
onds of  data  must  be  available  after  the  detected  phase  arrival  time.  Where 
no  detection  occurs,  at  least  25.  6 seconds  of  data  must  be  available  beginning 
at  the  expected  phase  arrival  time. 

In  the  event  that  fitting  a compressed  edit  is  not  possible,  the 
event  in  question  is  flagged  unavailable,  and  processing  skips  to  the  next 
event. 

After  edit  compression  is  completed,  certain  locations  in  the 
event  header  pertaining  to  edit  status  (i.  e.  , the  edit  start  time,  edit  length, 
and  detected  phase  arrival  time)  are  updated.  This  is  the  last  step  of  the  edit 
compression. 


2.  Measurement  of  Broadband  Complexity,  Mean  Phase  Standard 
Deviation  and  Maximum  Mean  Frequency 

Acceleration  data  are  obtained  by  differentiating  the  displace- 
ment trace.  At  this  point,  a broadband  complexity  discriminant  is  obtained 
for  teleseismic  events  (delta  > 20  degrees).  The  complexity  is  computed  as 
the  ratio  of  the  area  under  the  complex  envelope  in  the  first  10  seconds  im- 
mediately after  signal  onset  to  the  area  under  the  envelope  in  the  first  5 sec- 
onds. This  measure  is  used  instead  of  the  standard  definition  given  in  Sec- 
tion II  in  order  to  minimize  path  scattering  effects  which  appear  for  longer 
time  windows,  and  to  reduce  the  bias  arising  from  the  inclusion  of  large 
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noise  contributions,  especially  for  smaller  events.  Two  additional  discrimi- 
nants, maximum  mean  frequency  and  mean  phase  standard  deviation  are  also 
measured  from  the  acceleration  trace.  These  discriminants  are  discussed 
in  some  detail  by  Unger  (1978a),  and  in  Section  n. 


Both  the  maximum  mean  frequency  and  the  mean  phase  stand- 
ard deviation  are  obtained  from  the  quadratic  regression  estimate  of  the  con- 
tinuous instantaneous  phase.  An  arbitrary  time  signal,  x(t),  may  be  written 


in  terms  of  its  envelope,  Ex(t),  ant*  its  unwrapped  instantaneous  phase,  4*x(t). 


as  (Whalen,  1971): 


x(t)  = Ex(t)  cos  <*>x( t) 


(HI-1) 


where 


Ex(t)  = (x2(t)  + x2(t))" 


4>  (t)  = tan  1 x(t)/x(t)  . 


Here,  x(t)  is  the  Hilbert  transform  of  x(t)  , given  by 


x(t)  = PV 


f *LiL 

J *(t-r) 


dr 


where 


'PV'  denotes  the  Cauchy  principal  value. 


To  obtain  the  maximum  mean  frequency  and  mean  phase  stand- 
ard deviation,  a quadratic  regression  is  performed  to  fit  the  instantaneous 
phase  with  a quadratic  polynomial  within  a window  of  given  length  (Unger, 
1978a): 


*x(t.  r ) = aQ(t)  + a^tjr  + a2(t)r‘ 


where 
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t is  the  window  start  time 

t is  the  time  relative  to  the  start  of  the  window 

t)  is  the  estimated  instantaneous  phase  at  time 

t for  the  window  starting  at  time  t 
*0(t).a  (t).a  (t)  are  the  regression  coefficients  for  the  window 
beginning  at  time  t . 

In  the  regression  process,  $x(t)  is  fit  to  the  unwrapped  instantaneous  phase, 
^x(t,  r),  under  a minimum  mean-squared  error  criterion.  The  window  is 
then  stepped  ahead  one  sample  and  a new  regression  is  performed.  A 4 sec- 
ond window  is  used  in  this  study.  The  phase  standard  deviation  for  each  win- 
dow is  given  by: 

L 

Vt}  = 1 - #(t,r.)j2/(L-3)^  , 

i=l  L J 

where  L is  the  length  of  the  window  in  samples.  The  time  variant  mean 
phase  is  defined  as  the  regressed  phase,  evaluated  at  the  center  of  each  win- 
dow (t=  L/2): 


^x(t  + L/2)  = $x(t,  L/2)  = aQ(t)  + ax(t)L/2  + a2(t)L2/4.  (m-2) 

The  time  variant  mean  frequency  is  then  found  by  differentiating  equation 
(DI-2)  with  respect  to  r,  and  evaluating  at  r=  L/2: 

?(t  + L/2)  = £aj(t)  + La2(t)J  jin. 

The  broadband  complexity,  mean  frequency,  and  mean  phase 
standard  deviatiqn  are  only  obtained  if  the  broadband  signal  detector  is  trig- 
gered. Where  no  detection  is  declared,  the  measurement  of  these  quantities 
is  bypassed  and  a null  operation  flag  is  set  for  that  particular  event- station 
combination. 
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Short- Period  Magnitude  Measurements 


Several  short-period  magnitude  estimates  are  obtained  in  the 
course  of  short-period  event  processing.  For  teleseismic  events  (A  > 20°), 
a broadband  mfe  and  several  variable  frequency  magnitudes  (VFM)  are  cal- 
culated. These  are  scaled  using  the  Veith  and  Clawson  (1972)  P- factors.  In 
addition,  where  data  are  available  for  regional  events,  crustal  Pn,  Sn,  and 
Lg  magnitudes  are  calculated  according  to  the  following  formulas: 


m 


= 3.  82  + 2.  00  log  A + log  (A/T)  (Evernden,  1967) 

mSn  = 3-  79  + !•  8°  log  A + log  (A/T)  (Fitch  et  al.  , 1 978) ) (III-  3) 

mLg  = 3.  30  + I.  66  logA  + log  (A/T)  (Nuttli,  1973). 


The  logic  flow  for  the  broadband  amplitude  and  period  measure- 
ment is  shown  in  Figure  III- 6.  The  first  step  consists  of  several  quality  con- 
trol checks  concerning  the  availability  of  a short-period  < dit  and  the  consist- 
ency of  the  desired  measurements  with  available  data.  For  example,  checks 
are  made  to  ensure  that  regional  phase  measurements  are  not  performed  on 
teleseismic  date,  and  that  the  supplied  data  falls,  at  least  partially,  within 
the  following  velocity  windows: 


Pn 

Sn 

Lg 

Teleseismic  P 


8.  3 km/sec  - 5.  0 km/ sec 
5.  0 km/sec  -3.8  km/sec 
3.  8 km/sec  - 2.  4 km/sec 
15.  0 km/sec  - 8.  0 km/sec. 


If  these  checks  are  not  passed,  control  returns  to  the  supervisory  routine. 


If  the  checks  are  passed,  a 5 percent  cosine  taper  is  applied 
to  the  beginning  and  end  of  the  data  series.  The  series  is  then  broadband 
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FIGURE  III- 6 

BROADBAND  SHORT-PERIOD  (A/T)  MEASUREMENT  LOGIC  FLOW 
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filtered  between  0.  3 Hz  and  one-half  the  Nyquist  frequency.  Broadband  filter- 
ing is  performed  by  an  eighth  order  Butterworth  filter. 

After  the  data  trace  has  been  filtered,  the  signal  envelope  and 
instantaneous  frequency  are  generated.  The  signal  envelope  is  given  in  equa- 
tion (HI- 1 ),  where  the  Hilbert  transform  is  calculated  in  the  frequency  domain 
using  the  relation: 


x(f)  = -i  sgnf  x(f ) , 


where  sgnf  is  the  'signum'  operator  defined  by 


sgnf  = 


-1  f < 0 . 


The  instantaneous  frequency  is  obtained  by  evaluating  the  derivative  of  the 
instantaneous  phase  in  closed  form  from  equation  (HI- 3): 


where 


x(t)  = E^(t)  cos  4>x(t) 


<#>  (t)  = tan  1 x(t)/x(t) 


Since  the  instantaneous  frequency  f(t)  is  related  to  the  instantaneous  phase 


<*>x(  t)  = 2 it  f(t)dt, 
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Substituting  the  expression  for  4>x(t)  given  in  equation  (III -1)  and  differentiat- 
ing with  respect  to  t : 


3T 


^ (tan*1  x(t)/x(t)) 


x(t)^-  - *(t) 


<**(*) 

dt 


E*<‘> 


= 27Tf(t). 


(Ill- 4) 


The  expression  in  equation  (HI-4)  has  the  advantage  that  it  may  be  calculated 
without  first  unwrapping  the  instantaneous  phase.  A difference  operator  is 
used  to  approximate  dx(t)/dt  and  dx(t)/dt  . 

After  the  envelope  and  instantaneous  phase  have  been  obtained 
the  envelope  is  searched  for  the  maximum  amplitude  of  signal  or  noise. 
Whether  it  is  signal  or  noise  is  determined  by  whether  the  broadband  detector 
war  triggered  or  not.  If  the  detector  was  not  triggered,  the  entire  envelope 
is  searched  for  a noise  peak,  log1  nA^^/T  is  calculated  where  A^  is 
the  maximum  of  the  noise  envelope  and  T is  the  corresponding  period  deter- 
mined from  the  instantaneous  frequency,  and  this  value  is  stored  in  the  event 
header.  Where  the  detector  was  triggered  the  envelope  is  searched  in  a 10 
second  gate,  beginning  2 seconds  before  the  detected  signal  onset  time.  Again, 

log  ^A  / T is  calculated  for  the  signal  peak,  and  this  value  is  stored  in  the 

°10  max 

event  header.  The  log  A /T  values  are  later  scaled  to  teleseismic  or  re- 

® max 

gional  crustal  phase  magnitudes  using  the  formulas  given  in  equation  (HI-3). 

In  addition  to  the  broadband  magnitudes,  short -period  variable 
frequency  magnitudes  (VFM)  are  also  generated.  The  VFM  measurement 
logic  flow  is  illustrated  in  Figure  HI- 7. 
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SHORT- PERIOD  VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  MEASUREMENT  LOGIC  FLOW 


The  VFM  measurement  phase  begins  with  a series  of  quality 
control  checks  similar  to  those  used  for  the  broadband  magnitude  measure- 
ments. In  the  same  way,  if  any  of  these  checks  are  not  passed,  a null  mea- 
surement flag  is  set  and  control  is  returned  to  the  supervisor. 

The  short-period  VFM  ratios  are  calculated  in  eight  frequency 
bands  for  data  sampled  at  20  Hz,  where  the  center  frequencies  are  0.  32  Hz, 

0.  50  Hz,  0.  79  Hz,  1.  26  Hz,  2.  00  Hz,  3.  16  Hz,  5.  01  Hz,  and  7.  94  Hz.  For 
10  Hz  data,  only  six  frequency  bands  are  used  where  the  maximum  center 
frequency  is  3.  16  Hz.  The  S narrowband  filter  (described  in  Appendix  A) 
is  used,  scaled  to  constant  Q over  all  bandwidths.  The  bandwidth  to  center 
frequency  ratio  is  set  to  0.452,  which  allows  the  -6  dB  points  of  the  filter 
to  just  overlap,  ensuring  minimum  leakage  between  adjacent  frequency  bands. 
From  calibration  studies  it  was  also  found  necessary  to  scale  the  filter  output 
for  short- period  data  by  dividing  the  resultant  amplitudes  by  the  area  under 
the  filter,  given  by 


* E Hj‘ 

j=i 


where  H^  are  the  frequency  domain  filter  weights  and  Af  is  the  frequency 
domain  sampling  interval. 

Once  the  narrowband  trace  has  been  generated,  the  signal 
envelope  and  instantaneous  frequency  are  generated  in  the  same  manner  as 
for  the  broadband  measuremmts.  The  next  step  is  a test  to  determine 
whether  or  not  the  broadband  signal  detector  triggered  for  the  event  under 
consideration. 

Where  the  detector  did  not  trigger,  the  entire  envelope  is 
searched  for  a noise  peak  with  a corresponding  period  within  the  bandlimits 
of  the  filter.  If  no  noise  peak  is  found  within  those  bandlimits  (as  may  be  the 
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case,  due  to  leakage  of  strong  peaks  from  adjacent  frequency  bands)  a null 
measurement  flag  is  set  and  processing  continues  with  the  next  frequency 
band. 
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If  a noise  peak  with  a corresponding  period  within  the  filter 

bandlimits  is  found,  log(A  ) is  calculated  where  A is  the  maximum 

max  max 

noise  amplitude.  This  value  is  then  stored  in  the  event  header,  flagged  as  a 
noise  measurement. 

Where  the  broadband  signal  detector  did  trigger,  the  envelope 
is  searched  for  a maximum  in  a 10  second  window  starting  2 seconds  prior 
to  the  signal  onset  time.  If  no  signal  peak  with  a corresponding  period  with- 
in the  filter  bandpass  is  found,  the  frequency  band  is  declared  a noise  band 
and  a noise  measurement  is  attempted  as  described  previously.  If  a signal 
peak  is  found,  then  a detection  ratio  is  calculated  for  that  frequency  band  by 
taking  the  ratio  of  the  peak  signal  envelope  to  the  average  noise  envelope  in 
whatever  leading  noise  window  is  available.  The  log^^A^^)  for  the  signal 
peak  is  then  calculated  and  stored  in  the  event  header,  together  with  the  de- 
tection ratio,  and  processing  proceeds  to  the  next  frequency  band. 

As  in  the  case  of  the  broadband  magnitude  measurements,  the 
log,  „ (A  ) values  for  the  VFM  are  scaled  to  magnitude  using  the  formulas 
given  in  equation  (III- 3).  An  additional  threshold  test  is  employed  for  each 
frequency  band,  where  the  signal  envelope  peak  to  average  noise  envelope 
ratio  is  tested.  The  purpose  of  this  test  is  to  avoid  indiscriminant  signal 
declarations  in  frequency  bands  which  are  dominated  by  noise.  In  order  to 
maintain  signal  status,  the  signal  envelope  peak  must  exceed  the  average 
noise  envelope  by  at  least  10  dB.  Otherwise,  the  measurement  in  that  band 
is  treated  as  noise. 

In  addition  to  measuring  magnitudes  in  the  narrow  frequency 
bands,  a narrowband  complexity  measurement  is  also  made  in  each  frequency 
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band.  This  measurement  is  only  made  when  the  broadband  detector  is  trig- 
gered, and  is  carried  out  in  a manner  similar  to  that  used  in  obtaining  the 
broadband  complexity. 

C.  LONG- PERIOD  EVENT  MEASUREMENTS 

The  long-period  discriminant  processing  segment  of  TISSPROG 
consists  of  two  parts.  The  first  part  involves  removal  of  the  instrument  re- 
sponse, which  is  identical  to  the  procedure  used  to  correct  for  the  instrument 
response  for  short-period  data.  The  second  part  consists  of  the  timing  of  the 
long-period  phases  and  the  measurement  of  the  long-period  VFM.  This  is 
illustrated  schematically  in  Figure  III- 8. 

Because  no  long-period  signal  detector  is  available  to  operate 
reliably  on  long-period  data  without  a priori  specified  dispersion  data  and 
because  the  events  in  the  test  data  base  may  exhibit  vastly  different  dispersion 
characteristics  as  a result  of  the  numerous  different  travel  paths  to  various 
stations,  an  alternate  method  of  long- period  signal  timing  and  detection  is 
needed.  To  deal  with  this  problem,  broad- region  dispersion  models  com- 
posed by  Unger  (1978b)  are  used.  These  models  are  based  upon  group  velocity 
data  obtained  by  Sun  (1977).  The  group  velocity  data  consist  of  recordings  at 
a variety  of  stations  of  eastern  Kazakh  (EKZ)  and  Nevada  Test  Site  (NTS)  pre- 
sumed nuclear  explosions.  For  these  events,  it  is  possible  to  define  an  ex- 
ponentially modulated  group  velocity  window  for  any  frequency  of  interest,  as 
shown  in  Figure  HI- 9 for  NTS  events  and  Figurfe  III- 10  for  EKZ  events.  Al- 
though the  group  velocity  windows  are  of  necessity  quite  broad,  they  do  prove 
to  be  useful  in  timing  the  long-period  Love  and  Rayleigh  phases. 

Because  most  events  in  the  test  data  base  originate  in  the  Asian 
continent,  and  to  minimize  the  problem  of  running  off  the  end6  of  the  supplied 
long-period  edits,  the  eastern  Kazakh  timing  function  is  used  in  this  6tudy. 
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For  Love  wave  timing  at  a given  frequency  f , the  high  and  low  velocity 
limits  are  given  by 


V (f)  = 3.  0 + 0.  75  e'40(f'°' °2)  + 0.  4 

H 

VT  (f)  = 3.  0 + 0.  75  e"40(f‘0- 02)  - 0.  4 


Similarly,  the  Rayleigh  wave  velocity  limits  are  given  by 


V (f)  = 2.  75  + 0.  75  e'4°(f'°'02)  + 0.  4 

ri 


VL(f)  = 2.  75  + 0.  75  e"40(f'0,  02)  - 0.  4 . 


In  generating  the  long-period  VFM  measurements,  this  timing 
algorithm  is  used  to  re-edit  the  input  data  at  each  of  the  seven  frequencies  of 
interest.  In  a manner  similar  to  the  short-period  VFM  measurements,  the 
data  traces  for  the  long-period  vertical  and  long-period  transverse  components 
are  narrowband  filtered  about  center  frequencies  of  0.  02  Hz,  0.  025  Hz, 

0.  0316  Hz,  0.  04  Hz,  0.  05  Hz,  0.  0631  Hz,  and  0.  0794  Hz.  As  in  the  short- 
period  VFM  measurement  phase,  the  narrowband  filter  described  in  Appendix 
A is  used.  Again,  the  filter  has  constant  Q over  all  frequencies,  with  a 
bandwidth  to  center  frequency  ratio  of  0.  444.  This  allows  the  filters  to  over- 
lap at  slightly  less  than  the  -3  dB  points. 


The  first  step  in  the  long- period  VFM  measurement  involves 
running  the  long-period  timers  to  re-edit  the  data  at  each  filter  center  fre- 
quency of  interest.  The  Rayleigh  wave  timer  is  run  on  the  vertical  component. 
Next,  checks  are  made  to  ensure  that  the  new  edit  indeed  fits  within  the  sup- 
plied edit.  If  the  first  1. 5 cycles  of  data,  as  delimited  by  the  timing  function, 
are  not  available,  a null  measurement  flag  is  set  and  another  component  or 


another  frequency  band  is  processed.  If  the  first  1.  5 cycles  of  data  or  more 
are  available,  the  trace  is  re-edited  with  up  to  400  seconds  of  leading  noise 
included,  and  then  narrowband  filtered,  and  the  signal  envelope  and  instan- 
taneous frequency  are  computed  in  a maimer  identical  to  that  used  in  the 
short-period  VFM  measurement  phase. 

The  next  step  consists  of  a search  for  the  peak  in  the  signal 
envelope  (not  including  the  leading  noise  portion)  with  a corresponding  period 
within  j the  bandwidth  either  side  of  the  center  frequency.  Where  no  signal 
peak  is  found,  a null  measurement  flag  is  set  and  processing  proceeds  to  the 
next  component  or  next  frequency  band.  Where  a signal  peak  is  found,  a de- 
tection ratio  of  the  signal  envelope  peak  to  the  average  noise  envelope  is  form- 
ed. This  ratio  is  then  tested  against  a 10  dB  threshold.  Where  the  detection 

ratio  exceeds  the  threshold,  a signal  is  declared,  and  log, „(A  ) is  com- 

° *10  max 

puted,  where  A is  the  peak  signal  amplitude.  This  value  is  placed  in 

max 

the  event  header  and  is  later  scaled  to  magnitude. 

For  edits  wholly  contained  within  the  supplied  signal  edits, 

where  the  detection  ratio  fails  to  exceed  the  threshold,  the  measurement  is 

declared  noise  and  log1  n(Aryiay)  is  computed  and  placed  in  the  event  header, 

where  A is  now  the  peak  noise  amplitude.  In  the  case  of  short  edits 
max 

(i.  e.  , edits  not  wholly  contained  within  the  supplied  signal  edit),  a null  mea- 
surement flag  is  set  and  processing  proceeds  to  the  next  component  or  next 
frequency  band.  The  reason  this  is  done  instead  of  declaring  a noise  mea- 
surement is  to  avoid  the  possibility  of  declaring' a noise  measurement  for  a 
given  event  where  a large  signal  would  have  been  present  had  the  entire  edit 
been  obtained.  This  eliminates  the  possibility  of  assigning  a low  noise  esti- 
mate to  what  would  have  been  a large  detected  signal  had  the  edit  been  com- 
plete. 

In  addition  to  the  individual  VFM  magnitudes,  an  analogue  to 
a broadband  M is  developed  by  taking  the  maximum  of  the  long-period 
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vertical  VFM  measurements.  This  is  then  incorporated  into  the  standard 


Mg-  m^  discriminant. 


The  formulas  used  in  scaling  long-period  log(A/T)  to  magni- 
tude are  given  by  Nuttli  and  Kim  (1975)  for  Eurasian  earthquakes  and  explo- 


sions; 


Ms  = log10(A/T)  + 0.89  log10A  + 1.34 


A < 10 


Ms  = log10(A/T)  + !•  07  1°8iqA  + !•  16  10°<A<25° 


Mg  = log10(A/T)  + 1.66  log1()A  + 0.  34 


A > 25 


D.  EXAMPLES  OF  SHORT- PERIOD  AND  LONG- PERIOD  EVENT 


MEASUREMENT 


In  this  subsection,  two  examples  are  presented  which  illustrate 
the  typical  processing  of  short-period  and  long-period  records.  The  event 
processed  is  an  earthquake  which  occurred  in  eastern  Russia: 


Origin  time 


Latitude 


Longitude 


m (NEIS) 

D 


Depth  (NEIS) 


1 November  1977 


03:54:24 


55.  3°N 


130.  8°E 


33  km. 


The  short-period  data  used  here  were  recorded  at  Bluff, 
Alaska  (BFAK).  The  long-period  data  were  recorded  at  the  ASRO  site  at 
Kabul,  Afghanistan  (KAAO). 
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1.  Example  of  Short- Period  Event  Measurement 

The  computer  output  for  the  short-period  event  measurement 
example  is  shown  in  Figure  III- 11.  The  first  step  consists  of  reading  the  in- 
put seismogram  card  to  identify  the  event  to  be  processed  (Figure  III- 11a). 
The  data  are  fetched  in  segments  of  128  points,  decimated  if  desired,  and 
the  power  within  each  segment  is  computed.  After  the  data  trace  has  been 
fetched,  the  broadband  short-period  detector  is  applied  to  the  data  as  de- 
scribed in  Subsection  III-B.  Next,  the  system  response  i6  removed  as  de- 
scribed in  Subsection  III-B  (Figure  Ul-llb).  The  edit  is  then  compressed 
about  the  detected  (broadband)  phase  arrival  if  the  detector  has  triggered, 
or  about  the  expected  phase  arrival  if  the  detector  has  not  triggered  (Figure 
III-llc).  A plot  showing  the  original  edit  for  this  example,  together  with  the 
detected  phase  arrival  time  and  the  51. 2 seconc  compressed  edit  is  shown 
in  Figure  ID-12. 

After  edit  compression  has  been  completed,  acceleration  and 
displacement  traces  are  generated.  Measures  of  the  broadband  complexity, 
mean  phase  standard  deviation,  and  mean  instantaneous  frequency  are  then 
obtained  from  the  acceleration  trace  (Figure  ID- lid),  as  described  in  Sub- 
section III-B.  The  broadband  log1Q(A/T)  measurement  is  next  (Figure 
Dl-lle)  followed  by  the  VFM  log1()A  and  narrowband  complexity  measure- 
ments (Figure  ID-1  If ).  The  broadband  logjg(A/T),  the  VFM  log^A,  and 
narrowband  complexity  are  all  measured  from  the  displacement  trace.  Fig- 
ure ni-13  shows  the  narrowband  filtered  traces’used  in  the  VFM  log^A 
measurement  process,  together  with  plots  of  the  corresponding  envelopes 
and  instantaneous  frequencies.  The  final  step  in  the  short-period  event  pro- 
cessing is  the  organization  of  the  measurements  and  other  event  information 
into  the  event  header,  which,  together  with  the  compressed  data  edit,  is  then 
written  onto  the  output  tape  (Figure  ID-llg). 
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FIGURE  HI- 11a 

SHORT- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
SHORT- PERIOD  DATA  EDIT  AND  BROADBAND  SIGNAL  DETECTION 
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FIGURE  El- lib,  c 

SHORT- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
(b)  INSTRUMENT  RESPONSE  REMOVAL,  AND 
(c)  SHORT- PERIOD  EDIT  COMPRESSION 
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FIGURE  HI- lid,  e 

SHORT- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
(d)  BROADBAND  COMPLEXITY,  MEAN  INSTANTANEOUS 
FREQUENCY,  MEAN  PHASE  STANDARD  DEVIATION, 

AND  (e)  BROADBAND  LOG.  n(A/T)  MEASUREMENT 
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FIGURE  Hl-llf 

SHORT- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  LOG  A MEASUREMENT 
AND  NARROWBAND  COMPLEXITIES 
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FIGURE  in-llg 

SHORT- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
INFORMATION  CONTAINED  IN  SHORT- PERIOD  EVENT  HEADER 
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2.  Example  of  Long- Period  Event  Measurement 

The  output  for  the  long-period  event  measurement  example  is 
shown  in  Figure  III-14.  As  in  the  case  of  the  short-period  example,  the  first 
step  in  the  long-period  event  measurement  process  consists  of  reading  the 
input  seismogram  card  to  identify  the  event  to  be  processed  (Figure  III-14a). 
The  three  components  of  data  are  fetched  in  segments  of  128  points,  the  data 
are  rotated  to  the  vertical,  transverse,  and  radial  configuration,  and  auto- 
mated checks  are  made  to  determine  whether  the  data  contain  evidence  of 
clipping  or  uncorrectable  spikes.  Next,  the  instrument  response  is  removed 
(Figure  ni-14b),  and  measurements  of  log^A  are  made  on  the  vertical  and 
transverse  components  for  long-period  VFM  computation  (Figure  III- 14c). 

The  narrowband  filtered  long-period  traces  and  the  corresponding  envelopes 
and  instantaneous  frequencies  are  shown  in  Figure  III- 15.  As  in  the  short- 
period  case,  the  final  step  in  the  long-period  event  processing  is  the  organiza 
tion  of  the  long-period  measurements  and  other  event  information  into  the 
event  header,  which,  together  with  the  three  components  of  data,  is  then 
written  onto  the  output  tape  (Figure  III-14d). 
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FIGURE  ID- 14a 

LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING  LONG- PERIOD 
THREE  COMPONENT  DATA  EDIT  AND  AUTOMATED  QUALITY  CHECKS 
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FIGURE  HI -14b 

LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 

LONG- 

PERIOD  SYSTEM  RESPONSE  REMOVAL 
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LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING  LONG- PERIOD 
VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  LOGvqA  MEASUREMENTS 
ON  VERTICAL  (I)  AND  TRANSVERSE  (2)  COMPONENTS 

(PAGE  1 OF  4) 


HI- 44 


tni  rugiiiKi  mot  e.oiu  to  miu  h hi  tbiocxtv  units 

;;i  fiiniin«Fiuil»,iii(iilTiiVMa,> 

TNt  IM  SOXT  STITTS  0 SROIII  INTO  TN|  MINIM  MIT 

TNvse  its  }i«  s icc nos  or  foxti  if  tns  most  or  rot  soxt 

TNT  NM  fl*XT  STMTS  0 SBCONDS  INTO  Ml  OtOTSCN  MTT 


• ITS  MS  100  • IS.  SO  SIC 
08T0CTX2S  nr  XO  t«  CO  • 20.  4$ 


tic  rirooiNCT  •mat  o.ojh  ro  o.om  «i  ms  rsiocrrr  units 

TNt  kll  WIT  STMTS  0 SBCONDS  INTO  TNS  MTNMT  SOXT 
TN8II  ITT  ISO  SBCONDS  OT  NOXST  IT  TNK  MONT  Of  TNS  TOtT 

tni  mo  trr*  s-m-s  o sbcooos  into  tnb  obotbcn  boit 


2.02  ro  •.  »2  nt ntc 


••NCTTt  TNXS  XS  I SNOOT  BOXT** 

TNB  MI00MCT  NINOB  0.0000  TO  0.0000  NS  NIS  TBLOCXTT  UNITS 

HI 

TNT  NNO  OCXT  STINTS  0 SBCONDS  INTO  TNB  MXNIBT  BBXT 
TBBBT  MB  112  SBCONDS  Cf  NOXSI  IT  TNI  MONT  ON  TNB  BBXT 
TBB  MO  BDXT  STITTS  0 SBCOSOS  INTO  MB  OBOTBCB  BBXT 


loo#  ro  o.onm  as  bib  rvioem  units 

vn.,?r,R>wm 

• inni  mi  nt  mini  iu? 

ir  .out  it  tii  nm  of  m tttr 

• ttcoioa  mo  in  norm  tm 


LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING  LONG- PERIOD 
VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  LOG10A  MEASUREMENTS 
ON  VERTICAL  (1 ) AND  TRANSVERSE  (2)  COMPONENTS 

(PAGE  2 OF  4) 


^ — I ill  I -iii  — • -"g 

- _ _ . . " 


c. 


confoiim*  • • 


••moth  tmis  is  i irr  fbj?** 
tit  rrigmiT  fakof  o.osoo  ro  0.0100  ns  bas 

INF  CC»?»!JJf#  Ft  IT  MSS  Ml  lOUlfO  FOISTS 
tM  M|«t  HTlMfSt  FOUR  Of  T«0  IS  VII 

IS*  MTS  BBIT  STSF-S  0 SBCOBBS  IBTO  Til  FRIMAFT 
mil  IFF  m SFCOFOS  OF  MOIST  I?  *M  FtOMT  Of  TM  MBIT 
T»»  M«  FBIT  SHIFTS  0 StCOKBS  IMTO  TM*  GBOTBCM  id? 


CCITI  LUIFS  I.IT  *9  I.BI  FM/SSC 


0. w^QOO  0 ••  9. 012200 


'tw.tr 


it  timf-  m.o 


10-.19  -f  If  • I.OI 

VCIfF  If?  *k  V • H.M 


• ITM  FFRXCB  • 11.0%  SSC 
OTTfCTlOB  RATIO  IB  BB  • IS. OR 


:m***  »mio(  • ib.oo 

• OF  FT  I IM  DIM  At*  MR 


co«f:f»«*  • • 2 


?*•  msmcf  •a»of  o.osoo  to  o.osoo  m*  mas  vrlocitv  limits 

?MT  T9*FF*SSFB  F T IT  NAS  RAT  MOMIRRO  FOlRTS 
lM  B|f*  utiMIt  FOtFF  Of  TBO  IS  MI 

7RB  Bit  f TIT  STARTS  0 SFC0B03  IBTO  TM*  FRIMABT  (BIT 
Z M IF*  AF1  HI  S*CCMPS  OF  MOISF  AT  TM*  F*OBT  0*  TMI  *01* 

IMF  AC*  *riT  STARTS  J SBCOR0S  I*TO  *Mf  0*01 CC»  FBIT 


I.SM  TO  J.MR  * M/S  FC 


9.919000  ■ • 


D. 022200 


p\H  s*:  to  ••• 

1*0. 9*  FT*T'*I>«  22.2  S*C  IT  TUN*  M«.0 


10919  ••  if  i*fii-rrr  • 2.1* 

MIM  AVI  *MV  • 12. TA 


•ITM  FIRIOM  • 22.  20  SIC 
BCTCrrlO"  RATIO  IB  BB  • 22.  •• 


CfBTtA  F*FTOC  • 20.00 

• OF  FTS  IB  MOISI  At*  A MV  • 


re-F^RfM*  A • 


••1CT*:  T MI 5 IS  > • M9FT  IBIT** 

TMt  FFFOTFBCf  A A KOI  0.0111  TO  0.0AI1  R|  BAS  RFLOCITT  11FITS 

?w snnnu.!tiMi,or'iiDiiiivRiMtm 

Til  111  TBIT  STARTS  0 SBCOBBS  IBTO  TBF  F*I HAS V (BIT 
TMt FT  AFT  112  SFCOVPS  OF  MOISF  AT  TMt  FlOBT  OF  TMt  IBIT 
TMI  BIN  KBIT  ST# FTS  B SBCOBBS  IBTO  T0B  8B0TBCB  IBIT 


1.0  TO  l.M  «*/SBC 


9. D 199 JS  ?•  0.011100  ■ • 

Ml  Ft  OT  S FT  II  ft  A 
MAI*  11. 11  FFFtOr 


0.02B0IA 
IT.  T f FC 


IT  T IMF"  SI A.O 


19810  r>F  IF  A1FIITBPI  • I.Tf  VITA  FBRIOB  • 17. 7 

MOISI  Iff  tMV  • 0.07  BFTBCTfOR  RATIO  IB  Bl  • 

•••••••••••••••••••••••••••••••••MM 


SBC 

10.  17 


CBRTtF  FSIX3P  • IS. 01 

• OF  FTS  !■  BOISB  AH  BBT  • 11 


eONFON*BT  • • 2 


•MB  FRBOOBRC1  BARO*  0.9011  TO  0.0A11  IS  BAB  VBLOCITT  LIMITS 

TMI  IBI  BBIT  STARTS  0 SBCORBS  IBTO  TBt  FIIHABT  IBIT 
TBBBI  Atl  21«  SFCOB6S  OF  BOISB  IT  TMB  FlOBT  OF  TBt  B0IT 
TBB  111  BBIT  BTAFTI  B SBCOBBS  IBTO  TBB  OBOIBCM  IBIT 


l.M  10  1.H  IB/SBC 


11M  0.010MB  c*  0.01)190  B*  9.02I0H 

AM  FlfT  SI?  12  Ml 

FFAR"  li.n  FAFIOB"  17.  | SBC  AT  TIMB*  SBR.0 


19810  Zt  IF  A RFITTSOB  ■ 1.11  «IFB  FBBIOB  • 17.11  BBC 

■Ollt  At*. tMV  • 11.01  FTTBCTIOB  BITIO  IB  BB  • 1A.B0 

•••••••••••••••••••a 

FIGURE  HI- 14c 


• of  m ia  Mm  in  t>.  • *» 


L ONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING  LONG- PERIOD 
VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  LOGjqA  MEASUREMENTS 
ON  VERTICAL  (1 ) AND  TRANSVERSE  (2)  COMPONENTS 

(PAGE  3 OF  4) 


III-46 


L 


+9?  * *• 


mi  oircf  o.n»«  tc  o.otsr  hi  mis  volocxtt  units 

tut  eoNForsfr**  rr*T  mas  soi  rormid  ooirts 
TUX  RFF*  NTSflrS?  *0**0  or  r*n  is  512 

r*»!t  n»m  ii  seconds  tiro  mf  nmn  box? 

I"lOt  AM  «00  SVCCS**?  nr  <IOXS*  IT  TRt  F*0*T  Of  TM  *01? 

to*  ***  f*r?  stitts  m seconds  into  rur  codticn  «?*?? 


2. IT  ?0  1.*7  FF/SFr 


l3.«T  ’MIod*  11.*  src  IT  TIM* 


1*510  or  IF  >4Fi:T0OI  • 
KOXFF  AVF  ri>V  > A.I) 


* rr h Fcirco  • n.««  see 

DtTFCTIOR  RATIO  X N 0*  • 1*. 


7?**iF  OMTOf  ■ 12.5* 


or  fts  u roxsb  avi  r« 


•**0TF:  *015  :»  I SNOOT  I0XT** 

•■*  •FTO'IFNCT  MAC*  0. 07**  TO  0.07M  Ht  MS  fftOCITV  UNITS  2.* 

TM  CO*F»F**FO  FDIT  NAS  512  IOI1IIO  OOIRTS 
T 1*  OTI?  NTCNIST  FO*FF  Of  T*0  TS  512 

THF  MM  r.rx?  STAFFS  0 SBCORDS  XOTO  Tit  MINIM  10XT 
TNFFF  IFF  11*  SFCOODS  0?  NOXSF  IT  Ti«  FOOOT  CF  TNF  *DXT 

r«»  *r*  * pit  stiots  o stcovos  rro  tm  acortCR  nt* 

*2*  ;•  0.07**00  «•  0.01525* 

•••  OIOS  S*T  1ft  ••• 

FBI*-  *2.70  Ft F IOt*  1*.1  tie  IT  TI*t«  1*0.0 

tooio  or  tr  iNFiTTODi  • t.ti  ritb  obotoo  • j*. 15  src 

■OIBB  **«  *■*  • 5.00  DBTITTXO*  RATIO  IB  BO  • 17.0* 


2. *7  TO  1.7?  RN/SFT 


CPRTBI  riBIOP  • 12.5* 


FIGURE  DI-14C 

LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING  LONG- PERIOD 
VARIABLE  FREQUENCY  MAGNITUDE  (VFM)  LOGjqA  MEASUREMENTS 
ON  VERTICAL  (1 ) AND  TRANSVERSE  (2)  COMPONENTS 

(PAGE  4 OF  4) 


III- 47 


I 


n srr'Aositi*  mum 

J C3HWIIITJ  TATA  ll(CT«t  511  M»«« 

r»»ir  uhmumoi  amit  iuh  sio 

Dl?l  MWH  IWmi  }MI  H?f| 

1 1DITBC  SXTRf 

Hem  pait:t:op  rm  nomxi 
io«»:i  loimit  ssproc 

Mtl  Oi  XUTATIOP  t M* 

:ui3S'.r,u;»iirKfi?«:  i 

iMuriio  mu  1.0  loin 

Nil  ?IM|  15IJ*  |?rc»  77J05  ITPOAVI  -01- 

DATA  IIMTH  101*. 0 SACONrS 
SS9IAI  STASIS  AT  POT *T  10.  -1*0 

SOOPCi  TUtl  1«0*«  (ATCI  77105  (VN0AVI  -01' 

COPPXDtPCS  OP  fOnict  TXM  (TOP  CO0P| ! 

501PCP  AT  LATXT10T  $5. JO  I.  iO»CXT»»DI  110.00  I. 
A<T**ci«c  TO  OPT  . 5 • 0.5  0 • 0.0  TI  fll 


TX  1STX5AMD  SX0PAL-TO-KOXSS  PATIO!  0.0  TO. 
ST®.  0»f.  OP  ASS1D0AL  TUP:  STAT1 

»«o»i  • nrrpi  sossta  sot-proms  • i 

p-tm*  ?*•!?!  s-tise  lo-rwr  n- 


ItlWTN  (STATION  TO  SOOPCPU  01.7*  (PBINAPT  NAN  OX'KTXO')! 
0I9TANCB  (9*ATXON  TO  SCftPCS)  i Si.Sl  (®PC.|  -OP- 

STATION  ILPTATTON!  0.0  BN. 

StXSOONBTtO  TIPS:  SBCCSPlNft  ?TPt! 

T SC  TOPIC  CUSS  COOl!  AS 

LOO  APPlTTVOtSl  TSATTCAl  T*ANS tffPP 


LONG- PERIOD  EVENT  MEASUREMENT  EXAMPLE  SHOWING 
LONG -PERIOD  EVENT  HEADER  INFORMATION 


B W = 0.  009  Hz 

Vertical  Transverse 


v 

eu 

► 

c 

w 


w 


u 

e 

v 

3 

D* 

« 

ki 

A 

to 

3 

0 
V 

n 

1 

c 

(4 

•4J 

» 

a 


H 

XXX 

(M  o oo 

N N « 

o o o 

6 o o 


* 

I-  fi 
« 

& 

V 

A 
to 
3 

► i 

I 

rt 


ac  a sc 


fVJ 

CM 

O 


O 

N 

O 


00 

o 

d 


in 


W 

OS 

X 

A 


w to 
0 g W 

Bg5e 


S W r 
5 os  £ 
ZDO 


2 

< 


JDO- 
2 W X ^ r: 

K < CO 

5 2 S 50  W 

Q ^ £ z 

2 o<  ^ u 

. x ,n  ^ ^ 

CO 

y 

A 

o 

X 

w 


w 


8 » s § 

ShqS 

2 2 p 

2 W < to 
£ 2 n < 
>0  .wp 

A 

q 5 a y 

O g u < O 

OS  u 


1=1  < 


OS  O w -r* 

a‘*?89§§ 

2 2 s § “t 
a § s s o 

10  J 2 O ” 


Q 
W 
0S 

w 

H 

x 

M 

A 

Q 
2 

ffl  2 


a 

Q 

W 

to 

X 

to 

a 

o 


m P , 

2 a y 

gssg 

sis 

2 W n 

< H u 
, 2 W 

^ W 2 

< u s 

O ^ H 
S 06 
H w 


W N Q, 

!*  W W 

05  Q _ 

« w § 

g A 


06  h 

> 


X 

*-4 


III- 49 


Vertical 


5 ^ 

H {*}  | < 

5 ^ 0 C 
w I ^ 

§<JES 

5 2<qw 


w 3 H Q S 
W ^ S 2 g 
‘ft  5 w < W 

O £ g CQ  < _ 

w>g~§r 

> Q 5 h ^ h 

2 O S U < O 

W 5 1 1 °* 

w U [ft  »-H  ^ 

3‘Sga8g 

< o P5  2 J n] 

2 2 M W « ~ 

S3  j g a 
Q K < W 2 
w O « rK  W 
K U<  H fc  a 
W o n K ~ 

5|1  S5 

EDjZW 
w >4  w 

§«SUP 

«upg 

*g|5 

ga^E 

P5  W § 

5 05  O 


SECTION  IV 

MULTIVARIATE  DISCRIMINATION 


The  purpose  of  this  section  is  to  discuss  a statistical  process 
called  Multivariate  Discrimination, and  techniques  for  applying  that  process  to 
a large  Event  Discriminant  Data  Base  (EDDB),  in  an  effort  to  solve  the  event 
identification  problem.  The  section  begins  with  an  overview  presenting  the 
methodology  of  applying  multivariate  discrimination  to  the  event  identification 
problem.  It  continues  with  a detailed  description  of  the  EDDB  and  its  creation 
and  maintenance  as  an  integral  part  of  the  Extended  Interactive  Seismic  Pro- 
cessing System  (ISPSE)  on  the  PDP- 15/50  computer.  In  addition  the  analysis, 
implementation,  and  execution  of  a discrimination  algorithm  to  yield  event 
identifications  are  discussed.  The  algorithm  is  implemented  by  means  of  the 
Interactive  Seismic  Programming  Language  (ISPL)  which  facilitates  access  to 
the  EDDB  (Shaub  et  al.  , 1977).  This  section  concludes  with  a series  of  fig- 
ures illustrating  a typical  session  on  the  Event  Identification  System. 


OVERVIEW 


The  Event  Identification  System  was  created  for  the  purpose  of 
automatically  generating  signal  waveform  edits  and  measurements  to  provide 
a data  base  for  calculating  event  discriminants  and  for  classifying  the  type  of 
event  source  indicated  by  the  discriminant  measurements.  Since  a ground 
rule  of  this  task  was  to  classify  the  event  without  prior  knowledge  of  the  dis- 
criminant characteristics  of  a particular  source  region,  that  is  without  any 
training  data,  the  identification  system  was  designed  to  adaptively  isolate 
abnormal  events  and  build  up  a class  of  such  similar  events  as  other  events 
of  similar  discriminant  pattern  become  available. 


Past  research  has  indicated  that  most  discriminants  are  highly 
influenced  by  the  tectonic  region  of  the  source.  By  reciprocity,  it  might  be 
presumed  they  are  also  influenced  by  the  tectonic  region  of  the  receiver  since 
the  source  and  receiver  tectonics  convey  information  on  the  crustal  and  upper 
mantle  structure  traversed  by  the  raypaths.  For  example,  both  Sax  (1976) 
and  Unger  (1978a)  showed  that  the  pattern  of  short-period  discriminants  were 
completely  different  for  presumed  explosion  events  from  eastern  Kazakh, 
which  is  a tectonically  active  mountainous  region,  and  the  Nevada  Test  Site, 
which  is  located  in  a rift  zone.  To  take  care  of  this  situation  in  the  event 
classification  system,  all  signals  of  an  event  are  characterized  in  the  header 
by  a class  code  to  indicate  the  tectonic  classification  of  the  source  and  re- 
ceiver. This  can  be  used  either  to  derive  corrections  to  the  signal  measure- 
ments or  to  derive  separate  tectonic  class  population  means  and  standard  de- 
viations of  the  discriminants  derived  from  signals. 

A provision  is  made  in  the  system  to  visually  inspect  wave- 
forms for  the  purpose  of  checking  the  automatic  algorithms  used  to  window 
and  detect  signals.  As  a result  of  the  visual  inspection,  the  data  of  an  event 
can  be  prepared  again  to  assure  appropriate  selection  of  windows  and  detec- 
tion threshold.  Also,  based  on  selected  manual  runs  of  the  signal  measure- 
ments, the  signal  measurement  files  and  detection  status  of  the  measure- 
ments can  be  edited.  This  quality  control  process  is  considered  an  import- 
ant feature  of  the  system. 

The  generation  of  event  discriminants  from  magnitude  measure- 
ments of  signals  observed  at  a set  of  stations  is  done  by  applying  the  method 
of  Ringdal  (1974).  The  method  yields  an  unbiased  maximum  likelihood  esti- 
mate of  the  mean  magnitude  from  measurements  of  detected  and  non-detected 
signals.  Other  signal  measurements  are  processed  to  obtain  an  event  repre- 
sentation of  the  station  discriminants  derived  from  signal  measurements.  The 
resultant  event  discriminants  are  deposited  as  new  data  in  the  event 


discriminant  files.  In  case  of  changes  in  the  methodology  for  computing  a 
discriminant,  old  discriminant  data  can  be  retrieved,  corrected,  and  re- 
deposited. The  function  process  for  classifying  the  event  source  is  an  adap- 
tive part  of  the  event  identification  system.  The  classification  procedure  is 
started  by  determining  the  average  discriminant  vector  and  its  standard  devia- 
tion for  a set  of  'normal'  shallow  earthquake  events.  Each  unknown  event  is 
compared  to  the  average  earthquake  discriminant  vector.  If  the  event  is 
anomalous,  it  is  compared  to  other  anomalous  events  which  are  sub-divided 
into  sets  based  on  a measure  of  similarity.  If  it  does  not  fit  into  any  of  the 
prior  determined  sets  of  anomalous  events,  it  is  held  as  a singular  event 
until  enough  events  are  obtained  to  set  up  a new  class  of  anomalous  events. 

The  discrimination  capability  is  assessed  by  correlating  the  groupings  with 
known  or  presumed  earthquakes. 

Figure  IV- 1 illustrates  the  event  identification  system  as  four 
function  processes  operating  on  common  permanent  disk  files  of  seismic  data. 
The  files  are  accessed  by  any  user  of  the  system.  The  classification  status 
of  the  events  in  the  files  evolve  as  more  data  are  entered  into  the  system  and 
as  the  discriminant  data  files  are  subjected  to  more  classification  processing. 
The  classification  is  objective  given  prescribed  threshold  criteria  for  declar- 
ing an  event  anomalous  and  for  grouping  anomalous  events.  The  function 
which  prepares  discriminant  data  by  automatically  editing  events  and  detect- 
ing signals  is  performed  on  the  IBM  360/44  computer.  The  functions  which 
generate  event  discriminants  and  classify  the  event  sources  are  performed 
on  the  PDP- 15/50  computer. 

B.  THE  EVENT  DISCRIMINANT  DATA  BASE 

The  structure,  creation,  and  maintenance  of  the  Event  Dis- 
criminant Data  Base  (EDDB)  is  discussed  in  this  subsection.  The  EDDB  is 
resident  on  the  ten  megaword  disk  drive  of  the  PDP- 15/50  computer  and  is 
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ACTIONS 

(1)  Deposit  waveform  edits  on  disk  for  visual  inspection. 

(2)  Spot  check  the  detection  status  of  signals.  Revise  the  detection  status,  if 
necessary. 

(3)  Spot  check  the  signal  edit  windows.  Re-prepare  data,  if  necessary. 

(4)  Deposit  automatically  determined  signal  measurements  into  permanent 
seismic  measurement  files. 

(5)  Reduce  signal  measurements  to  a set  of  unbiased  estimates  of  the  event's 
discriminant  vector. 

(6)  Deposit  estimates  of  the  event's  discriminants  into  the  event  discriminant 
files. 

(7)  Classify  the  event  type,  if  possible,  otherwise  update  its  process  status 
and  try  again  when  more  event  data  are  available. 

FIGURE  IV -1 

EVENT  IDENTIFICATION  SYSTEM 


supported  by  the  Extended  Interactive  Seismic  Processing  System  (ISPSE) 
which  executes  on  that  machine  (Shaub  et  al. , 1977).  The  principal  items 
comprising  the  EDDB  are  vectors  of  discriminant  measurements  obtained  by 
preprocessing  waveform  data  from  various  recording  stations. 

1.  Description  of  the  EDDB 

The  general  organization  of  the  Event  Discriminant  Data  Base 
(EDDB)  is  shown  in  Figure  IV-2.  Note  that  it  consists  of  four  file  types, 
namely: 

• An  Event  Directory 

• Raw  Event  Discriminant  Data  Files 

• Averaged  Event  Discriminant  Class  Files 

• A Class  Directory  File. 

These  files  are  physically  compatible  in  that  they  are  easily  accessed  via  the 
DREAD  and  DWRIT  statements  featured  by  the  Interactive  Seismic  Program- 
ming Language  (ISPL).  (Subsection  IV-C-5,  The  Extended  Interactive  Seismic 
Processing  System  (ISPSE),  Shaub  et  al.  , 1977.  ) 

EDDB  structure  is  based  primarily  on  three  attributes:  event 
sequence  number,  station  designation,  and  tectonic  class  code.  When  an 
event  is  recorded  in  the  Event  Directory  it  i6  assigned  a sequence  number. 

For  each  event  so  recorded,  at  the  most  fifty  discriminant  vector  entries, 
corresponding  to  fifty  reporting  stations,  reside  in  one  or  two  Raw  Event  Dis- 
criminant Data  Files.  The  Event  Directory  groups  these  entries  by  pointing 
to  Raw  Event  Discriminant  Data  File  locations  holding  discriminant  vectors 
for  stations  belonging  to  the  same  tectonic  class. 

Currently,  there  are  three  tectonic  classes  in  the  EDDB;  'A'  = 
active,  'I'=inactive,  and  'R'  = rift.  Hence,  the  tectonic  class  code  for  an  event 
with  epicenter  in  an  'A'  zone  which  is  recorded  at  a station  in  an  'I'  zone  is 
'Al'.  Codes  derived  in  this  way  are  also  called  event- station  class  by  class 


codes  and  are  the  names  given  to  the  Averaged  Event  Discriminant  Class 
Files.  These  files  contain  average  discriminant  vector  entries,  one  per 
event,  computed  as  an  unbiased  mean  (Ringdal,  1975)  of  the  discriminant 
vectors  for  all  reporting  stations  belonging  to  the  same  tectonic  class.  There 
fore,  if  an  event  is  reported  by  stations  spanning  all  three  tectonic  classes,  it 
will  appear  in  chree  Averaged  Event  Discriminant  Class  Files. 


To  unify  the  above  concepts,  consider  the  remaining  EDDB 
file  type,  the  Class  Directory.  The  Class  Directory  is  partitioned  by  the  nine 
tectonic  class  codes.  Within  each  partition,  the  corresponding  Averaged 
Event  Discriminant  Class  File  name  is  recorded,  together  with  a vector  of 
Event  Directory  sequence  numbers.  This  vector  specifies  which  events  be- 
long to  that  file  and  their  locations  within  the  file. 


With  an  overview  of  the  EDDB  structure  in  mind,  detailed  de 
scriptions  of  the  four  file  types  noted  will  now  be  presented. 


The  structure  of  the  Event  Directory  is  shown  in  Figure  IV- 3 
A representative  entry  to  the  Event  Directory  is  comprised  of  attributes 
which  may  be  defined  as  follows: 


Event  Sequence  Number  - An  integer  J,  1<J<199,  assigned 
to  an  event  when  it  is  recorded  in  the  Event  Directory.  For 
example,  the  fifth  event  recorded  is  assigned  the  event  se- 
quence number  5. 


Event  Designation  Number  - The  unique  four  digit  number 
appearing  as  the  last  four  characters  of  the  TI  alphanumeric 
event  designation. 


Event  Class  Code  - The  tectonic  class  ('A',  'I', 
the  zone  in  which  the  event's  epicenter  is  located 


• Station  Class  Code  and  Status  Vector  - A fifty  character  string 
where  each  character  represents  the  tectonic  class  ('A1,  'I',  or 
'R')  of  the  associated  reporting  station.  The  character  0 in- 
dicates that  the  station  did  not  report  the  event. 

• Event  Station  Class  by  Class  Code  Pointer  Vector  - A vector 
whose  elements  point  to  Raw  Event  Discriminant  Data  File  loca- 
tions holding  discriminant  vectors  that  have  the  same  event- 
station  class  by  class  code.  Elements  are  partitioned  by  at  the 
most  three  different  codes. 

In  addition  to  the  above  information,  the  Event  Directory  also 
contains  a header  record  specifying  the  total  number  of  events  in  the  directory, 
the  number  of  Raw  Event  Discriminant  Data  Files  utilized  for  discriminant 
vector  entries,  and  the  largest  Event  Sequence  Number  appearing  in  the  Class 
Directory.  These  parameters  are  necessary  for  EDDB  maintenance  functions 
to  be  discussed  in  this  subsection. 

b.  Raw  Event  Discriminant  Data  Files 

All  the  Raw  Event  Discriminant  Data  Files  have  the  same 
structure  as  shown  in  Figure  IV-4.  A representative  Raw  Event  Discriminant 
Data  File  entry  is  comprised  of  attributes  which  may  be  defined  as  follows: 

• Discriminant  Vector  - Array  of  discriminant  measurements 
obtained  by  preprocessing  waveform  data  from  various  record- 
ing stations.  The  array  is  partitioned  in  four  categories  de- 
pending on  the  type  of  data  from  which  the  discriminants  are  ' 
derived,  namely,  long-period  signals/noise,  short-period 
regional  signal/noise,  short-period  teleseismic  signals /noise, 
and  short-period  signals/noise. 

• Discriminant  Status  Vector  - Five  element  array  of  ones  and 
zeros  indicating  that  discriminants  derived  from  various  data 
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types  are  present  or  absent  in  the  discriminant  vector.  The 
data  types  are  defined  as  follows;  long-period  signal,  long- 
period  noise,  short-period  signal,  short-period  noise,  and 
regional  teleseismic. 


Event-Station  Class  by  Class  Code  - Tectonic  class  code  de- 
rived by  combining  the  event  class  code  with  the  tectonic  class 
of  the  reporting  station,  e.  g.  , 'AI1. 


Event  Name  - The  TI  alphanumeric  event  designation 


Station  Name  and  Number  - The  TI  alphanumeric  station 
designation  together  with  its  assigned  integer  number  K, 

1 < K < 30. 


Event  Origin  Time  - Source  time  in  'YY  DDD  HH  MM  SS1 


Event  Latitude  and  Longitude  - Source  latitude  and  longitude 
in  degrees,  N(+),  S(-),  E(+),  and  W(-). 


Event  Depth  - Source  depth  in  km 


Event  magnitude 


Averaged  Event  Discriminant  Class  Files 


All  the  Averaged  Event  Discriminant  Class  Files  have  the  same 
structure  as  shown  in  Figure  IV-5.  A representative  Averaged  Event  Dis- 
criminant Class  File  entry  is  comprised  of  attributes  which  may  be  defined 
as  follows: 


Average  Discriminant  Vector  - The  unbiased  mean  (Ringdal, 
1975)  of  all  discriminant  vectors  having  the  same  event- station 
class  by  class  code  for  the  given  event. 


Cumulative  Discriminant  Status  Vector  - The  sum  of  all  dis 
criminant  status  vectors  having  the  same  event- station  class 
by  class  code  for  the  given  event. 


STRUCTURE  OF  AN  AVERAGED  EVENT  DISCRIMINANT  CLASS  FILE 


Averaged  Event  Discriminant  Class  File  Name  and  Number 
The  tectonic  class  code  = event- station  class  by  class  code 
associated  with  this  file,  e.  g.  , 'AI'. 


Event  Designation  Number  - The  unique  four  digit  number 
appearing  as  the  last  four  characters  of  the  TI  alphanumeric 
event  designation. 


Number  of  Raw  Discriminant  Vectors  Contributing  to  Average 
Discriminant  Vector  - Self  defining. 


Class  Directory 


The  structure  of  the  Class  Directory  is  shown  in  Figure  IV-6. 
Note  that  all  attributes  comprising  a representative  entry  to  the  Class  Direc 
tory  are  either  self-defining  or  have  been  previously  defined  in  connection 
with  other  EDDB  files. 


The  Event  Discriminant  Data  Base  (EDDB)  is  generated  from 
seven  track  Discriminant  Tapes  produced  by  the  seismic  data  preparation 
software  on  the  IBM  360/44  computer  (Schmidt,  1978).  A detailed  descrip- 
tion of  the  individual  discriminant  measurements  stored  on  these  tapes  and 
the  algorithms  from  which  the  measurements  are  derived  may  be  found  in 
Section  HI  of  this  report. 


The  off-line  PDP-15/50  program  DDBASE  reads  the  seven 
track  Discriminant  Tapes  and  initializes  the  EDDB  Event  Directory  and  Raw 
Event  Discriminant  Data  Files  discussed  in  Subsection  IV-B-1.  The  general 
flow  of  DDBASE  is  shown  in  Figure  IV- 7.  Note  that  the  DDBASE  initializes 
the  EDDB  one  event  at  a time,  hence  new  events  may  be  continually  added  to 
the  data  base  as  they  are  processed. 


Read  Event  Directory  Header 
to  Determine  EDDB 
Initialisation  Starting  Point 


FIGURE  IV- 7 

GENERAL  FLOW  OF  DDBASE 
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Maintenance  of  the  EDDB 


Subsequent  to  EDDB  initialization  by  means  of  DDBASE,  the 
data  base  is  maintained  by  user  functions  programmed  in  the  Interactive  Seis- 
mic Programming  Language  (ISPL).  Since  ISPL  is  an  integral  part  of  the 
Extended  Interactive  Seismic  Processing  System  (ISPSE),  the  user  functions 
take  advantage  of  the  system's  graphics  console  environment  to  facilitate  in- 
teractive retrieval  and  editing  of  information  contained  in  the  EDDB  and  inter- 
active construction  of  the  Averaged  Event  Discriminant  Class  Files  discussed 
in  Subsection  IV-B-1. 

a.  DLOGD 

The  user  function  DLOGD  enables  the  ISPSE  user  to  interac- 
tively list  (DELOG)  the  Event  Directory,  DELOG  and  display  discriminant 
measurements,  or  DELOG  and  edit  discriminant  measurements.  The  general 
flow  of  DLOGD  is  shown  in  Figure  IV- 8.  The  console  communications  from  a 
sample  DLOGD  execution  is  shown  in  Figure  IV- 9. 

b.  CLAVE 

The  user  function  CLAVE  enables  the  ISPSE  user  to  interac- 
tively construct  the  Averaged  Event  Discriminant  Class  Files  and  Class 
Directory  previously  discussed.  CLAVE  provides  the  option  of  combining 
tectonic  classes,  that  is,  to  average  the  discriminant  measurements  for  all 
stations  reporting  the  event  regardless  of  class  code.  CLAVE  also  permits 
the  user  to  specify  particular  stations  to  be  included  in  the  averaging  process. 
Note  that  all  magnitude  averages  are  computed  using  Ringdal's  method  for 
determining  unbiased  network  magnitude  estimates  (Ringdal,  1975).  The  gen- 
eral flow  of  CLAVE  is  shown  in  Figure  IV- 10.  The  console  communications 
from  a sample  CLAVE  execution  is  shown  in  Figure  IV-11. 
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SAMPLE  CLAVE  EXECUTION 


C.  APPLICATION  OF  DISCRIMINATION  OPERATORS  TO  THE  EVENT 
DISCRIMINANT  DATA  BASE  (EDDB) 

The  purpose  of  this  section  is  to  describe  an  event  identifica- 
tion algorithm  for  classifying  events  either  as  earthquakes  or  as  possible 
explosions,  and  to  provide  a description  of  the  software  designed  for  imple- 
mentation of  that  algorithm.  A multivariate  approach  is  taken  to  the  problem 
of  event  identification  in  this  study.  The  basis  for  the  event  identification 
algorithm  is  the  hypothesis  that  given  a statistically  representative  set  of 
normal  earthquakes,  it  is  possible  to  identify  anomalous  events  by  their 
deviation  from  the  mean  of  the  base  set.  To  this  end,  z- statistics  relative 
to  the  means  and  standard  deviations  of  the  discriminants  over  the  earthquake 
set  are  calculated,  and  correlation  terms  are  computed  between  events. 
Threshold  tests  are  then  used  to  classify  the  events  on  the  basis  of  these  cor- 
relation terms. 

The  Event  Discriminant  Data  Base  (EDDB)  management  sys- 
tem allows  the  analyst  to  look  at  the  edited  signal  phases  from  which  the  dis- 
criminant measurements  were  derived  for  a particular  event- station  pair. 
This  is  done  with  the  aid  of  the  system  functions  DPSCAN,  SELEV,  and 
FILTER  described  by  Ringdal  et  al.  (1975).  With  these  functions  it  is  pos- 
sible to  select  any  phase  desired  and  to  study  any  frequency  band  of  interest. 
It  is  also  possible  to  develop  user  functions  to  make  independent  measure- 
ments of  the  various  quantities  used  to  calculate  the  discriminants.  These 
capabilities,  together  with  the  (previously  mentioned)  EDDB  maintenance 
functions  discussed  in  Subsection  IV-B  provide  the  system  with  a large  de- 
gree of  flexibility.  Specifically,  an  analyst  may  easily  override  machine 
decisions  and/or  measurements  and  substitute  his  own.  Note  that  editing 
operations  via  DLOGD  should  be  done  before  performing  the  user  function 
CLAVE  to  calculate  the  average  unbiased  measurements  for  each  event. 


The  multivariate  analysis  user  functions  for  event  identifica- 
tion operate  on  the  average  class  files  which  have  been  generated  by  the  user 
function  CLAVE.  The  actual  discriminant  values  are  calculated  from  the 
average  event  measurements  and  employed  in  identifying  the  events.  The 
discriminants  and  their  positions  in  the  class  file  records  are  listed  in  Table 
IV- 1. 

The  general  flow  of  the  event  identification  algorithm  is  as 
follows.  The  first  step  is  to  select  a base  set  of  normal  earthquakes.  This 
base  set  is  used  to  identify  as  many  of  the  other  events  as  possible  as  belong- 
ing to  one  of  several  other  possible  sets.  These  sets  include  uncorrelated 
events,  possible  explosions,  and  several  types  of  presumed  explosions.  The 
characteristics  of  each  set  are  determined  by  its  initial  members.  In  other 
words,  the  various  sets  represent  the  properties  of  the  events  forming  the 
sets.  Classification  of  the  events  is  based  on  threshold  inner  products  (cor- 
relation terms)  of  vectors  for  each  event  containing  a z- statistic  for  each 
discriminant  value  computed  against  the  mean  and  standard  deviation  for  the 
discriminant  over  the  earthquake  set.  A z- statistic,  z,  is  defined  as 
z =(x- yd/cr  , where  x is  a given  observation  and  fi  and  a are  a related 
mean  and  standard  deviation,  respectively.  The  projections  for  a particular 
event  must  then  pass  or  fail  a threshold  test  a specified  number  of  times  be- 
! . fore  the  event  is  moved  to  other  than  the  uncorrelated  or  possible  explosion 

set.  As  the  presumed  explosion  sets  are  filled  attempts  are  made  to  sub- 
classify each  set.  Events  that  correlate  above  a preset  threshold  are  used 
to  generate  a mean  z-vector  for  the  set.  These  mean  z-vectors  are  used  in 
computing  later  correlation  terms.  Events  that  do  not  meet  this  criterion 
are  moved  into  the  next  explosion  set  where  further  subclassification  is  at- 
tempted when  enough  events  are  present.  In  this  manner,  an  event  sifts 
down  through  the  explosion  sets  until  it  is  either  classified  as  one  of  eight 
separate  types  or  it  passes  into  the  set  representing  all  other  types  of  ex- 
plosions. This  general  set  does  not  have  a mean  z-vector  and  subclassification 
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of  it  is  not  possible.  Thi6  procedure  is  intended  to  serve  only  as  a bootstrap 
for  the  event  identification  methods  used  by  Sax  (1976).  Ho  ver,  if  carried 
far  enough,  this  procedure  itself  can  be  used  for  event  identification.  The 
methods  used  by  Sax  involve  projection  of  vectors  of  earthquake- based  z- 
statistics  for  each  discriminant  against  a mean  z-vector  for  a given  type  of 
explosion.  The  mean  z-vectors  used  in  this  method  would  be  obtained  by  the 
method  described  above.  The  event  identification  bootstrapping  software  is 
described  below. 

The  user  function  DBOOT  provides  the  analyst  with  a system 
of  auxiliary  user  functions  which  operate  on  the  EDDB  with  the  goal  of  event 
identification.  The  general  flow  of  DBOOT  is  as  follows.  First,  the  class 
file  on  which  to  operate  is  selected  by  inputing  the  appropriate  class  number 
in  the  range  from  1 to  9.  If  the  selected  class  file  does  not  contain  any  events, 
DBOOT  is  terminated.  Following  this  the  various  parameters,  earthquake 
means,  and  standard  deviations  (SD),  and  the  partition  link  list  (defining  the 
event  sets)  are  restored  from  the  class  directory  to  their  last  state  (i.  e.  , as 
they  were  at  the  end  of  the  previous  identification  session  for  the  class  file). 
Various  user  functions  can  then  be  performed  at  the  analyst's  discretion  by 
making  selections  from  a descriptive  menu  of  functions. 

These  functions  are  described  in  detail  later  in  this  subsection. 
Briefly,  they  are  as  follows: 

• PARM  - enables  the  analyst  to  set  the  various  parameters 
used  throughout  the  other  functions. 

• SETEQ  - allows  a set  of  known  earthquakes  to  be  selected 
from  the  events  in  a class  file  to  serve  as  a base  for  classifi- 
cation. 

• SETU  - adds  new  events  to  the  uncorrelated  (U)  set. 
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• EQNRM  - calculates  the  means  and  standard  deviations  of 

the  discriminants  over  the  events  in  the  earthquake  set. 
Calculates,  ranks,  and  displays  a modified  norm  for  each 
event  vector  of  z- statistics  based  on  the  means  and  standard 
deviations. 

• PXGEN  - generates  the  possible  explosion  (PX)  set  from  the 

U set. 

• EXPCL  - attempts  to  subclassify  the  PX  set  and  to  generate 

a file  of  probable  explosions  (XI  set). 

• SUBEX  - attempts  to  subclassify  a particular  set  of  probable 

explosions  and  to  generate  a mean  z-vector  for  the  set. 

• PARTY  - displays  the  events  in  a given  set. 

Upon  termination  of  DBOOT  the  parameters,  earthquake  mean  and  standard 
deviation  vectors,  and  the  partition  link  list  are  saved  in  the  class  directory. 
Figure  IV- 12  illustrates  the  structure  and  flow  of  DBOOT. 

PARM  enables  the  analyst  to  set  the  various  parameters  used 
throughout  the  bootstrapping  procedure.  These  parameters  are  primarily 
comprised  of  various  test  thresholds  and  counter  limits.  A description  of 
each  parameter  and  its  previous  value  is  displayed  each  time  this  function 
is  performed.  Parameters  can  be  set  or  reset  at  any  time  during  classifica- 
tion; however,  in  general  it  is  recommended  that  the  parameters  be  kept  con- 
stant during  a given  classification  attempt.  Otherwise,  inconsistent  classifi- 
cations may  result.  Figure  IV- 13  illustrates  PARM. 

The  function  of  SETEQ  is  to  allow  a set  of  earthquakes  to  be 
selected  from  the  events  in  a given  class  file  to  act  as  a basis  for  further 
classification.  The  event  sequence  number  of  a 'known'  earthquake  is  input 
and  the  event  sequence  numbers  in  the  file  are  searched  for  this  number. 
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USER  FUNCTION  PARM 


If  it  is  found,  it  is  added  to  the  earthquake  set.  All  events  that  are  not  as- 
signed to  the  earthquake  set  are  assigned  to  the  uncorrelated  set.  After  all 
the  earthquakes  are  selected,the  earthquake  and  uncorrelated  sets  of  the  par- 
tition link  list  are  saved  and  all  other  sets  and  counters  are  zeroed.  This 
function  is  only  performed  in  order  to  reinitialize  the  classification  process. 
The  flow  of  SETEQ  is  illustrated  in  Figure  IV- 14. 

SETU  provides  the  capability  of  assigning  new  events  to  the 
uncorrelated  set  (U-file).  These  new  events  can  then  be  operated  on  by  the 
other  functions  available  in  DBOOT.  Figure  IV-15  illustrates  the  function 
SETU. 


EQNRM  generates  z- statistics  and  threshold  norms  for  the 
events  resident  in  the  earthquake  set.  First,  the  mean  and  standard  devia- 
tion are  computed  for  each  discriminant  across  the  earthquake  set.  The  z- 
statistics  are  then  computed  against  these  means  and  standard  deviations. 
The  ’norm'  of  the  event  is  the  square  root  of  a threshold  inner  product  of  the 
vector  of  z- statistics  for  the  event  with  itself.  The  threshold  inner  product, 

A-B,  A » (aj,  a2,  . . . , a ),  B ■ (b^,  b^,  ....  b^^),  may  be  written 
ND 

A-B  = v w.a.b.  where  ND  is  the  number  of  discriminant  z- statistics 
j=l  J J J 

and 


r 

w.  = < 

J ( 1 


if  la  ,|<  t or  |b.|<  t or  a.b.  < 0 
1 J1  1 J1  J J 


(IV- 1) 


otherwise 


where  r is  a selected  threshold  on  the  z- statistic.  This  type  of  correlation 
term  is  used  in  later  functions  as  a measure  of  the  similarity  of  two  events. 

It  emphasizes  the  major  similarities  between  events  and  allows  for  a certain 
degree  of  ’noise’  by  making  use  of  the  threshold  on  the  vector  components. 
After  the  norms  are  computed  for  all  of  the  events  they  are  ranked  from  high 
to  low  and  are  displayed.  These  norms  are  used  in  determining  the  test 
threshold  in  PXGEN  to  be  discussed  next.  Ideally,  the  threshold  should  be 
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significantly  inside  the  range  displayed  by  the  earthquakes  in  order  to  ensure 
that  explosions  are  not  erroneously  classified  as  earthquakes.  Finally,  the 
mean  and  standard  deviations  are  saved  before  returning  to  DBOOT  to  allow 
restoration  of  the  system  in  future  runs.  Figure  IV- 16  diagrams  the  flow  of 
EQNRM. 

PXGEN  generates  a set  of  possible  explosions  (PX  - file)  from 
the  U - file.  An  event  is  selected  from  a display  of  those  in  the  U file.  A 
z- statistic  is  computed  against  the  earthquake  mean  and  standard  deviation 
for  each  discriminant  associated  with  the  event.  A modified  norm  of  this 
z-vector  is  then  computed  by  taking  the  square  root  of  the  threshold  auto- 
inner  product.  This  norm  is  then  tested  against  a predetermined  threshold. 

On  the  basis  of  this  threshold  test,  the  event  is  classified  as  an  earthquake 
(£  threshold)  or  a possible  explosion  (>  threshold),  and  it  is  moved  to  the 
appropriate  set.  The  flow  of  PXGEN  is  diagrammed  in  Figure  IV-17. 

EXPCL  attempts  to  subclassify  the  possible  explosions  into  a 
'known'  explosion  group  providing  the  number  in  the  file  {PX  - file)  meets  or 
exceeds  a specified  minimum  set  in  PARM.  An  event  is  selected  from  a dis- 
play of  those  resident  in  the  PX  - file.  A threshold  inner  product,  or  cor- 
relation term,  is  then  computed  for  the  selected  event's  z-vector  against 
those  z-vectors  representing  all  the  other  events.  The  maximum  correlation 
term  is  then  chosen  and  tested  against  a preset  threshold.  If  the  event  pair 
associated  with  the  maximum  correlation  term  has  never  been  chosen  before 
or  if  the  test  result  is  different  than  before,  the  event  pair  is  considered 
valid.  Otherwise,  it  is  discarded  and  the  next  highest  correlation  term  is 
passed  through  the  same  tests.  This  procedure  is  repeated  until  either  a 
valid  pair  is  found  or  all  event  pairs  available  have  been  tested.  Once  a valid 
pair  has  been  found  the  results  are  stored  and  counters  are  changed  accord- 
ingly for  both  events  in  the  pair.  If  the  event  pair  has  passed  the  threshold 
» 

test,  that  is,  if  the  correlation  term  is  greater  than  the  threshold,  the  success 
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SD  = Standard  Deviation 


FIGURE  IV- 16 
USER  FUNCTION  EQNRM 
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counter  (C)  is  incremented  by  one  and  the  failure  counter  is  decremented  by 
some  weight  (set  in  PARM).  If  the  success  counter  meets  or  exceeds  a pre- 
viously specified  maximum,  the  selected  event  is  moved  to  the  first  explosion 
set  (XI  - file),  and  its  associated  counters  are  zeroed.  If  on  the  other  hand, 
the  event  pair  fails  the  threshold  test,  the  correlation  term  being  less  than 
or  equal  to  the  threshold,  the  failure  counter  (K)  is  incremented  by  one.  If 
the  failure  counter  meets  or  exceeds  the  specified  maximum  the  selected 
event  is  moved  to  the  earthquake  set  and  its  counter  is  zeroed.  Otherwise, 
the  event  is  moved  to  the  uncorrelated  set  and  re-execution  of  PXGEN  is  re- 
quired to  move  it  back  into  the  PX  - file  before  it  may  be  re-correlated.  The 
objective  of  decrementing  the  failure  counter  upon  passing  the  test  and  mov- 
ing an  event  to  the  uncorrelated  partition  upon  failing  is  to  help  prevent  an  ex- 
plosion from  being  misclassified  as  an  earthquake.  It  is  suggested  that  the 
minimum  number  of  events  necessary  to  attempt  subclassification  be  set 
significantly  higher  than  the  maxima  on  the  test  counters.  Otherwise,  erron- 
eous classifications  may  result  due  to  a predominance  of  marginal  earthquakes 
over  explosions.  Figure  IV-18  shows  the  logic  flow  of  EXPCL. 


SUBEX  is  analogous  to  EXPCL  in  that  it  attempts  to  subclassify 
a set.  In  the  case  of  SUBEX,  a particular  explosion  set  is  chosen  in  which  to 


i : 


attempt  subclassification.  As  before,  a specified  number  of  events  is  neces- 
sary before  subclassification  can  be  initiated.  After  this  condition  has  been 
met,  it  is  possible  to  select  an  event  from  a display  of  those  in  the  set.  How- 


ever, only  events  that  have  never  passed  the  threshold  test,  that  is,  events 
that  are  not  already  included  in  the  mean  z-vector  for  the  set,  can  be  select- 
ed. The  threshold  inner  product  or  correlation  term  is  again  computed 
between  the  selected  event's  z-vector  and  the  z-vectors  of  all  the  other  events 
If  an  event  has  previously  passed  the  threshold  test  either  as  a selected  event 
or  a matched  event  the  mean  z-vector  for  the  set  is  used  in  computing  the 
correlation  term.  The  maximum  correlation  term  is  then  subjected  to  a 
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FIGURE  IV- 18 
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threshold  test  to  determine  the  disposition  of  the  selected  event.  When  an 
event  pair  passes  the  test,their  z- statistic  vectors  are  used  to  update  the 
mean  z-vector  for  the  set  if  they  have  not  already  been  included.  As  pre- 
viously mentioned,  these  events  are  then  excluded  from  selection  by  the 
analyst.  When  a pair  fails  the  test,  the  failure  counter  (NN)  is  incremented 
by  one.  When  the  counter  reaches  or  exceeds  a preset  maximum  the  event 
is  moved  to  the  next  explosion  set  and  all  relevant  counters  for  the  event  are 
zeroed.  The  particular  characteristics  delineating  the  various  explosion 
types  is  determined  by  the  first  few  events  that  are  correlated.  Figure  IV-19 
details  the  flow  of  SUBEX. 

The  functions  that  have  been  described  above  are  the  backbone 
of  the  event  identification  system.  By  making  use  of  these  functions  the  anal- 
yst can  generate  a solid  foundation  for  classification  of  additional  events.  It 
should  be  noted  that  a great  deal  of  experimentation  might  be  necessary  be- 
fore the  various  parameters  involved  in  this  method  can  be  set  with  any  de- 
gree of  accuracy.  The  system  available  allows  the  analyst  a great  deal  of 
freedom  for  this  experimentation. 

D.  EXAMPLES  OF  DISCRIMINATION  OPERATOR  APPLICATION 

The  purpose  of  this  subsection  is  to  illustrate  the  execution  of 
the  event  identification  user  functions  described  in  the  previous  subsection. 
The  following  figures  and  annotations  typify  the  analyst- machine  interaction 
during  an  event  identification  session. 

Before  the  event  identification  session  can  be  initiated,  the 
ISPSE  system  must  be  executed.  The  procedure  for  doing  so  is  described 
by  Shaub  and  Black  (1977).  After  this  has  been  done,  the  event  identification 
system  is  initiated  by  performing  DBOOT.  This  routine  is  the  driver  for  the 
remaining  identification  functions.  Figure  IV- 20  shows  the  execution  of 
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FIGURE  IV- 20 
EXECUTION  OF  DBOOT 
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DBOOT.  Note  the  menu  describing  the  various  user  functions  that  can  be 
performed  from  DBOOT.  This  menu  is  displayed  after  every  execution  of 
one  of  the  subsidiary  functions. 

Figure  IV-21  illustrates  the  selection  of  'normal'  earthquakes 
by  using  SETEQ.  The  number  of  earthquakes  in  the  base  population  is  dis- 
played as  each  match  is  found.  SETEQ  can  be  aborted  without  affecting  pre- 
vious classifications  and  counters  by  inputting  a negative  number.  After  an 
earthquake  set  has  been  selected  the  function  is  exited  by  inputting  zero  as 
shown. 

Figure  IV- 22  illustrates  execution  of  PARTY.  This  function 
can  be  used  to  display  any  event  classification  partition  by  inputting  the  proper 
partition  number  as  seen  in  the  figure.  The  contents  are  listed  with  their  in- 
dex into  the  partition,  their  event  sequence  number,  and  their  event  designa- 
tion number.  The  event  sequence  number  is  the  numerical  order  of  the 
events  in  the  Event  Discriminant  Data  Base  (EDDB).  The  event  designation 
number  is  the  identifying  number  for  the  event  and  corresponds  to  the  event 
number  in  the  Area  of  Interest  event  listings.  In  this  example,  the  earth- 
quakes  selected  in  SETEQ  are  displayed. 


Figure  IV-23  illustrates  execution  of  SETU.  The  error  mes- 
sages shown  on  the  figure  result  from  trying  to  access  empty  partitions. 

This  occurs  when  determining  how  many  events  are  included  in  the  various 
groups.  As  seen  in  this  figure,  all  the  events  are  accounted  for  in  the  parti- 
tions. There  are  no  new  events  to  add  to  the  uncorrelated  file. 

Figure  IV- 24  shows  how  PARM  is  used  to  set  the  parameters 
used  throughout  the  identification  procedure.  Note  that  the  character  string 
'VALUE'  must  be  input  to  change  the  value  for  a particular  parameter  and 
the  string  'EXIT'  must  be  input  to  terminate  the  function. 
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Figure  IV- 25  illustrates  execution  of  EQNRM.  This  routine 
calculates  threshold  norms  for  vectors  of  discriminant  z- statistics  for  each 
earthquake.  The  method  is  described  in  the  previous  subsection.  The  sorted 
norms  seen  in  the  figure  can  be  used  in  selecting  a threshold  for  the  threshold 
test  in  PXGEN.  The  system  errors  also  present  in  the  figure  for  the  operator 
XSTAT  occur  when  the  discriminant  z- statistics  are  computed  for  an  event. 
They  indicate  one  or  more  unavailable  discriminants  in  one  or  more  events. 
The  validity  of  the  various  computations  is  not  affected.  Empty  discriminant 
components  are  due  to  the  unavailability  of  signal  windows  or  of  detected  sig- 
nals. The  calculator  of  the  norms  and  also  of  correlation  terms  (projections) 
is  designed  to  account  for  this  situation.  This  error  message  is  also  seen  in 
the  figures  relating  to  PXGEN,  EXPCL,  and  SUBEX  (Figures  IV- 26,  IV- 27, 
and  IV-28). 

Figure  IV- 26  illustrates  execution  of  PXGEN.  The  event  select 
ed  from  the  displayed  list  can  be  classified  as  either  an  anomalous  event  or  an 
earthquake  on  the  basis  of  its  vector  norm.  This  function  can  be  performed 
until  all  entries  in  the  uncorrelated  (U)-file  have  been  classified  as  either  a 
possible  explosion/anomalous  event  or  an  earthquake. 

Figure  IV-27  shows  how  EXPCL  is  used  to  subclassify  possible 
explosions  or  other  possibly  anomalous  events  into  other  classes  (earthquake 
or  explosion/anomalous  event).  As  explained  in  the  previous  subsection,  this 
subclassification  is  done  on  the  basis  of  threshold  tests  on  the  correlation 
terms  (projections  in  the  figure)  and  counters  on  the  numbers  of  successes 
and  failures.  In  this  example,  the  event  is  moved  into  the  first  explosion/ 
anomalous  event  partition.  At  present,  it  is  possible  to  move  all  events  ex- 
cept one  into  other  partitions  by  reiterating  this  function. 

Figure  IV-28  illustrates  execution  of  SUBEX  in  an  attempt  to 
subclassify  partition  number  4.  As  can  be  seen  in  the  figure,  only  events 
that  have  not  already  been  included  in  the  mean  z-vector  can  be  selected. 
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FIGURE  IV- 27 
EXECUTION  OF  EXPCL 
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The  selected  event  is  then  either  included  in  the  mean  z-vector,  as  in  this 
example,  or  moved  to  the  next  explosion/anomalous  event  partition.  As  ex- 
plained in  the  previous  subsection,  this  decision  is  also  based  on  the  basis  of 
threshold  tests  on  the  projections  and  counters.  The  matching  event  (i.  e.  , 
the  one  with  the  largest  projection  against  the  selected  event)  is  also  used  to 
update  the  mean  z-vector  if  it  has  not  been  included.  This  is  also  illustrated 
in  Figure  IV-  28. 

A possible  procedure  for  using  the  functions  illustrated  above 
to  identify  events  is  as  follows: 

• A base  set  of  earthquakes  is  selected  by  using  SETEQ. 

• Any  new  events  that  are  not  resident  in  the  partition  linklist 
are  added  to  the  uncorrelated  file  by  using  SETU. 

• The  threshold  on  the  calculation  of  the  norms  and  projections 
is  set  using  PARM. 

• EQNRM  is  used  to  generate  the  mean  and  standard  deviation 
vectors  for  the  earthquake  set  and  to  compute  and  display  the 
norms  of  the  z- statistic  vectors  for  the  earthquakes. 

• The  threshold  on  the  norms  (used  in  PXGEN)  is  set  using 
PARM.  It  is  also  possible  to  set  the  thresholds  and  counter 
limits  for  EXPCL  and  SUBEX  at  this  time  if  previous  exper- 
ience is  available.  Otherwise,  it  might  be  necessary  to  test 
run  several  of  the  classification  functions  and  reinitialize 
classification  with  SETEQ  in  order  to  gain  the  experience 
necessary  to  set  the  parameters. 

• PXGEN  is  used  to  classify  all  the  events  in  the  uncorrelated 
file  as  other  earthquakes  or  possible  explosions/anomalous 
events. 
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• The  possible  explosions  or  anomalous  events  are  further  class 
ified  using  EXPCL. 

• The  explosion/anomalous  event  partitions  are  then  evaluated 
using  SUBEX  to  determine  if  the  events  in  a particular  parti- 
tion are  similar. 

The  above  scenario  is  only  a rough  outline  of  the  various  pos- 
sible procedures  for  using  the  event  identification  functions.  It  should  also 
be  noted  that  the  system  described  and  illustrated  here  and  in  the  previous 
subsections  are  prototypes  and  are  subject  to  modifications. 


SECTION  V 
RESULTS 


An  event  identification  package  was  applied  to  a data  base  of 
thirty-five  events.  The  purpose  was  to  quality  check  all  components  of  the 
system,  to  obtain  experience  in  setting  thresholds  for  automatic  detection 
and  timing  of  edited  signals,  and  for  identification  of  normal  shallow  explo- 
sions and  other  event  types.  As  expected  as  a result  of  this  test,  we  identi- 
fied intermittant  program  malfunctions  and  uncovered  problems  in  designing 
signal  filters,  in  properly  timing  edits,  in  defining  independent  discriminants, 
in  properly  operating  the  classification  processor,  and  in  optimally  defining 
and  applying  multivariate  discriminants.  Our  results  indicated  the  feasibility 
of  our  goal  of  designing  and  operating  a seismic  event  identification  system 
which  automatically  edits  data,  computes  discriminants,  and  which  adaptively 
identifies  and  delimits  groupings  of  unusual  events  without  prior  knowledge 
of  such  events.  The  detailed  evaluation  of  these  results  are  given  in  the  re- 
mainder of  this  section. 

Table  V-l  shows  the  discriminant  calculations  used  in  this  pre- 
liminary analysis.  For  brevity  in  future  references,  each  discriminant  will 
be  referred  to  by  the  discriminant  number  shown  in  Table  V-l. 

The  first  five  discriminant  numbers  are  long-period  discrimi- 
nants. The  first  discriminant  indicates  the  average  center  band  magnitude 
of  Rayleigh  waves  between  16  seconds  and  32  seconds  compared  to  the  longer 
period  magnitude  between  40  seconds  and  50  seconds.  The  second  discrimi- 
nant measures  the  shorter  period  Rayleigh  wave  magnitude  at  1 3 seconds  com- 
pared to  the  Rayleigh  wave  center  band  magnitude.  The  third  discriminant 
measures  the  log-log  slope  of  the  Rayleigh  waves  between  13  seconds  and  16 
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seconds.  The  fourth  discriminant  measures  short-period  compressional 
wave  magnitude  relative  to  the  average  Rayleigh  wave  center  band  magni- 
tude representing  versus  reduced  to  a scalar  discriminant.  The  fifth 
discriminant  measures  the  Love  wave  magnitude  relative  to  the  Rayleigh  wave 
magnitude  at  a period  of  25  seconds. 

The  discriminant  numbers  six  through  twelve  are  short-period 
discriminants.  The  sixth  discriminant  measures  the  magnitude  at  0.  3 Hz 
relative  to  m^.  The  seventh  discriminant  measures  the  average  log-log  slope 
of  the  P wave  between  0.  3 Hz  and  0.  8 Hz;  the  eighth,  between  0.  5 Hz  and  0.  8 
Hz,  and  the  ninth,  between  0.  8 Hz  and  1.  3 Hz.  The  absence  of  frequencies 
higher  than  1.  3 Hz  was  caused  by  malfunction  of  filters  due  to  leakage,  which 
needs  to  be  corrected  before  valid  results  can  be  shown.  The  tenth  and  elev- 
enth discriminants  are  scalar  quantities  representing  the  relationship  between 
the  frequency  and  phase  standard  deviation  of  P waves.  The  twelfth  discrimi- 
nant is  a measure  of  complexity  of  P waves.  The  absence  of  discriminants 
based  on  measurements  of  regional  phases  was  due  to  the  sparsity  of  mea- 
surements available  in  the  data  base. 

Table  V-l,  part  2,  is  a listing  of  event  discriminant  measure- 
ments by  discriminant  number  and  event  numbers.  The  event  designation 
number  refers  to  the  event  numbers  in  Table  1-2  which  is  an  event  location 
list.  The  event  sequence  number  indicates  the  order  in  which  events  were 
entered  into  the  event  identification  disk  file.  The  event  list  in  Table  V-l 
is  grouped  to  show  the  initial  normal  shallow  earthquake  set  and  to  show  the 
initial  set  of  unknown  events  entered  into  the  identification  system. 

Figure  V-l  shows  a schematic  representation  of  a discriminant. 
The  population  on  the  left  represents  normal  shallow  earthquakes.  The  popula- 
tion on  the  right  represents  the  unusual  group  of  events  we  wish  to  correlate. 


the  events  in  the  order  of  the  measured  discriminant  values.  The  50%  level 
and  85%  level  of  the  ordered  statistics  are  obtained  in  the  sense  shown  on 
Figure  V-l.  The  deviation  parameter  | Diff  | is  obtained  as  the  absolute 
value  of  the  difference  between  the  85%  level  and  the  50%  level. 

Normal  shallow  earthquake  statistics  are  shown  for  the  initially 
selected  earthquake  training  set  on  Table  V-2.  For  each  discriminant  number, 
these  results  show  the  statistical  parameters  described  above  and  illustrated 
on  Figure  V-l.  The  multivariate  discriminant  statistics  are  shown  for  com- 
parison with  the  individual  discriminants.  How  the  multivariate  discriminant 
is  calculated  will  be  described  later.  Note  that  the  median  of  the  multivariate 
discriminant  is  a slightly  biased  statistic  for  the  ''ormal  earthquake  population. 
Table  V-3  shows  an  updated  statistical  evaluation  of  the  normal  earthquakes. 
The  update  was  obtained  after  adding  events  classified  as  earthquakes  to  the 
initial  earthquake  training  set.  Comparing  the  results  on  Table  V-3  to  those 
on  Table  V-2,  the  medians  appear  to  be  quite  stable  but  the  measure  of  devia- 
tion | Diff  | is  not  so  stable  possibly  indicating  that  more  data  should  have 
been  used  for  the  earthquake  training  set.  The  likelihood  ratio  shown  on 
Table  V-3  is  one  possible  measure  for  ranking  the  effectiveness  of  the  dis- 
criminants. Larger  numbers  indicate  more  effective  separation  of  the  non- 
normal group  from  normal  earthquakes.  It  is  computed  as  the  number  of 
events  in  the  correlated  set  above  the  85%  level  of  the  normal  earthquake 
population  divided  by  the  expected  number  if  the  group  were  normal  earth- 
quakes, that  is  0. 15  times  the  number  of  events  in  the  correlated  group.  If 
all  of  the  correlated  events  were  detected  by  the  discriminant,  the  likelihood 
value  would  be  6.  7.  That  was  the  case  for  the  multivariate  discriminant. 

Table  V-4  shows  the  discriminant  statistics  for  the  correlated  group.  These 
are  characterized  by  generally  larger  measures  of  deviation.  This  indicates 
that  misclassification  of  the  classified  group  as  normal  earthquakes  is  prob- 
ably a more  likely  error  than  misclassification  of  normal  earthquakes  as 
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TABLE  V- 2 

NORMAL  EARTHQUAKE  STATISTICS 


Discriminant 

Number 

Median 

85%  Level 

m 

1 

0.  36 

0.  57 

0.  21 

2 

-0.38 

-0.  26 

0.  12 

3 

-3.  97 

-3.  31 

0.  66 

4 

co 

o 

• 

1 

-0.  52 

0.  51 

5 

0.  05 

IT\ 

O 

O 

I 

0.  10 

6 

0.  63 

0.  19 

0.  44 

7 

-2.  11 

-0.  79 

1.  32 

8 

-2.  37 

-2.  09 

0.  28 

9 

-3.  01 

-2.  39 

0.  62 

10 

0.  01 

0.  05 

0.  04 

11 

1. 53 

1.  91 

0.  38 

12 

1. 44 

1.  32 

0.  12 

MVD* 

0.  95 

9.  23 

8.  28 

*MVD:  Multivariate  Discriminant 
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UPDATED  NORMAL  EARTHQUAKE  STATISTICS  WITH  THRESHOLD 
SET  AT  85%  LEVEL  OF  EARTHQUAKE  POPULATION 


Discriminant 

Number 

Median 

85%  Level 

IB 

Likelihood 

Ratio*** 

1 

0.42 

0.  59 

0.  17 

2.  6 

2 

1 

o 

CO 

UJ 

-0.  10 

0.  23 

3.  6 

3 

-3.  46 

-2.  39 

1. 07 

3.  1 

4 

-1. 08 

-0.  66 

0.42 

4.  6 

5 

-0.01 

-0.  07 

0.  06 

2.  1 

6 

-0.  82 

0.  21 

0.  61 

4.  1 

7 

-2.  36 

-0.  96 

1. 40 

5.  1 

8 

-2.47 

-1. 95 

0.  52 

5.  6 

9 

1 

UJ 

o 

OJ 

-2.  61 

0.  42 

3.  6 

10 

-0.  02 

0.  05 

0.  07 

6.  2 

11 

1.47 

1. 89 

0.42 

2.  1 

12 

1.43 

1. 26 

0.  17 

0.  6 

MVD* 

1.  22 

6.45 

5.  23 

6.  7** 

* MVD:  Multivariate  Discriminant. 

**  Maximum  attainable  value  for  this  likelihood  test  representing  no 
error. 

***  Likelihood  of  correct  identification  of  the  correlated  group. 
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TABLE  V-4 


CORRELATED  GROUP  STATISTICS  WITH  THRESHOLD 
SET  AT  85%  LEVEL  OF  EXPLOSION  POPULATION 


Discriminant 

Number 

Median 

85%  Level 

tea 

Likelihood 

Ratio*** 

1 

0.  58 

0.  26 

0.  32 

2.4 

2 

1 

o 

o 

■^3 

-0.  21 

0.  14 

4.  7 

3 

-2.43 

-4.40 

1. 97 

1.6 

4 

-0.  08 

o 

00 

o 

0.  72 

5.  5 

5 

-0.  02 

0.  21 

0.  22 

0.  4 

6 

-0.  12 

0.  55 

0.  67 

4.  3 

7 

-0.  11 

-2.  09 

1. 98 

3.  9 

8 

-1.42 

-1. 65 

0.  23 

5.9 

9 

-2.  14 

-3.  51 

1.  37 

2.  0 

10 

0.  39 

0.  25 

0.  14 

6.  7** 

11 

1. 64 

1. 47 

0.  17 

3.4 

12 

1.  36 

1.40 

0.  04 

3.  6 

MVD* 

56.  76 

37.  44 

19.  32 

6.  7** 

* MVD:  Multivariate  Discriminant. 

**  Maximum  attainable  value  for  this  likelihood  test  representing  no 

error.  i 

***  Likelihood  of  correct  identification  of  earthquakes. 
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members  of  the  correlated  group.  The  exception  is  discriminant  number 
ten  which  provided  complete  separation  of  the  normal  earthquake  below 
the  85%  level  of  the  correlated  group.  The  likelihood  ratio  for  correct  identi- 
fication of  earthquakes  was  calculated  as  the  number  of  normal  earthquakes 
beyond  the  85%  level  of  the  correlated  group  compared  to  the  expected  number 
of  correlated  events  beyond  their  85%  level,  i.  e.  , 0.  15  times  the  number  of 
events  in  the  correlated  group.  Note  that  the  multivariate  discriminant  in- 
dicated complete  separation  of  normal  earthquakes  below  the  85%  level  of  the 
correlated  group. 

Table  V-5  shows  event  discriminants  sorted  by  the  measured 
discriminant  values.  This  is  a convenient  form  for  representing  populations 
to  determine  the  medians  and  85%  levels.  Also,  the  event  designation  num- 
bers are  given  for  the  ranking  of  events  with  respect  to  each  discriminant. 

Based  on  the  results  of  the  event  classification  analysis,  six 
initially  unknown  events  are  correlated  with  normal  earthquakes.  These  were 
added  to  the  initial  earthquake  training  set  to  obtain  an  updated  earthquake 
training  set.  Event  designation  number  68  was  dropped  from  the  training  set 
due  to  an  edit  timing  error  caused  by  a programming  error.  The  event  iden- 
tification process  identified  fourteen  events  as  similar  members  of  a group 
designated  as  correlated  group  number  one.  In  that  group,  event  designation 
number  seven  was  dropped  due  to  a programming  error  in  the  timing  of  edits. 
Three  events  were  categorized  as  not  normal  earthquakes  and  as  not  members 
of  the  correlated  group  number  one.  Table  V-6  summarizes  these  results. 
Refer  back  to  Table  V-5  for  the  relative  ranking  of  events  in  these  groupings 
with  respect  to  particular  discriminants.  Refer  back  to  Table  V-2  through 
V-4  for  statistical  summaries. 

Table  V-7  shows  the  detectability  of  the  correlated  event  group 
number  one  and  three  uncorrelated  events,  which  are  at  this  point  consider- 
ed as  unknown  and  which  are  not  normal  shallow  depth  earthquakes.  The 
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TABLE  V-5 


EVENTS  SORTED  BY  MEASURED  DISCRIMINANT  VALUE 

(PAGE  1 OF  13) 


SORTED  EVENT  DISCRIMINANTS 


DISCRIMINANT 
EVENT  SEQ  4 

I. ^oogdo 

19.03030 

2i3.iV.J0OO 

II. 30000 

27.00000 
_ 29,03030 

6.  300003 
35.30000 
21.03000 

2.000000 
O.OOO00O 

33.00033 

26.00030 

9.000000 
33.00003 

22.00000 
10.00000 
5.030003 
13.000.10 
a. 003000 
25,0303(1 
15.033M0 
3.033003 

16.00030 
lO.P'igOn 
31.03100 

24. 00000 
•2.00000 

23.00000 
17.30330 
I'l.  0P0U.J 

34.00000 

20.00033 

7.000000 
32.00300 


*s  1.000003 
EVENT  DESGN  U 
36.03300 

59.00000 

14.00000 
3.330000 

21.00000 

17.00000 

73.00000 

53.03000 

61.00000 
48.03030 

19.03000 
16.00330 

77.00000 

20.30000 

36.00000 

62.00000 

1.003000 

66.20000 

34.00000 

56.00000 

151.0000 

37.00000 

55.00000 
39.00300 

50.00000 

61.00000 

149.0000 

7.000000 

143. 0000 

49.O003P 

35.00000 

47.30000 

22.03000 

10.00000 

30.00030 


D1SCR,  VALUE 
1.261485 
0.6961941 
0.6923602 
e. 8534517 
0,7 1788M l 
0.7027812 
0,6621035 
0.6238362 
0.6015817 
0,6004386 
0.5818741 
0,5760638 
0,5698195 
0.5604762 
0.5336183 
0,5171143 
0,5139959 
0.5119907 
0,5063442 
0.4516214 
0.4301729 
0.4246477 
0.3606163 
0.2931902 
0.2669533 
0.2585233 
0.2450056 
0.2186967 
0.2034087 
0. 1919603 
0.1635789 
0.1496212 
0.1121235 
0. 1051522 
0.665741OE-01 
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SORTED  EVENT  DISCRIMINANTS 


DISCRIMINANT 
EVENT  St'O  » 

as  2,003033 

EVENT  DESGN  » 

DISCR.  VALUF 

7.000030 

18.00000 

0.6141447 

2 9.DJI0.IW 

17.00000 

0.2291999 

27.00000 

21.00000 

0. 1520773 

I9.go.joa 

59,30000 

0.1419783 

2.O03O0O 

49.00300 

0,1000224 

11.00300 

4. 003003 

0.7S90890E-01 

33.00000 

36.00000 

-0.2436399E-01 

1.000330 

38,30000 

-0.2597964E-01 

10.00000 

t. 000000 

•0.4314971 E-0 1 

35,  JO0.10 

53.30000 

-0.6327271E-01 

30.00000 

16,00000 

-0.702444 3E-01 

21.00000 

61.00000 

-0.8165729E-01 

3.003030 

55.00000 

-0.1118910 

22.00000 

62.00000 

-0.1154736 

51.03003 

81.00000 

-0.1335688 

20.000-10 

14,30000 

-0.1445462 

4.000000 

56.00000 

-0.1457167 

0.000300 

19,00000 

-0.2102412 

9.000(130 

20.00000 

-0.2345200 

6.000000 

73.00000 

-0.2348187 

34.00000 

47.00000 

-0.2578360 

15.00000 

37,30000 

-0.2720613 

12.00000 

7.030000 

-0.3009868 

25.00000 

151.0000 

• -0.3326520 

13.30000 

34.30000 

-0.3596938 

24.00000 

149.0003 

•0.3614200 

5.000300 

68.00300 

-0.3996466 

14. 0-1030 

35.00300 

-0.4087738 

16..J3.J00 

39.03000 

-0.4637179 

2d. 00300 

22.00003 

-0.5044906 

17.00000 

49,00000 

-0.6350853 

18.0003(1 

53.00033 

-0.6417273 

32.33030 

30.00000 

-0.6679270 

26.00000 

77,00000 

-0.6705384 

23.33003 

143.0000 

-1.216998 
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EVENTS  SORTED  BY  MEASURED  DISCRIMINANT  VALUE 
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SORTED  EVENT  DISCRIMINANTS 


DISCRIMINANT 
EVENT  SEO  * 
29.00010 

11.000. 10 

34.01000 

17.00000 

35.03003 
1.1103000 

3.030000 

22.00003 

33.03000 

7.0000. 13 

27.00000 

31.03000 
2.000300 

21.00000 

4.0000P0 

33.00000 
6.001300 

32.00000 

17.00000 

14.00003 
5. 330303 
20.03030 

24.00000 
8.300030 
15. J00PU 

13.00000 

20.30003 
9.3 30000 

1.2.03000 

25.30000 

to. 00000 

10.30003 

18.03003 
26.00330 
23.103)0 


4s  3.000000 
EVENT  DESGN  * 

17.00000 

3.000000 
•47.00030 
59,0ni00 

53.00000 

38,03000 

55.00000 

62.00000 
16,30000 

18.00000 

21.00000 

81.01000 
•48.00000 

61,00000 

56.10000 

36.00000 

73.00000 

30.00000 
49,30030 
35.33000 

88.00000 

22.00000 

149.0000 

19.00000 

37.00000 

34.10000 

14.00000 

21.00000 

7.01.1000 

151.0000 

39.00000 

1.001000 

50 .00000 

77.00000 

143.0000 


OISCR,  VALUE 
0,8997595 
-0.3603995 
-0.8469236 
-1.357415 
-1,397402 
-1.885780 
-1,895353 
-2.094998 
-2.341238 
-2.377530 
-2.432505 
-2.507324 
-2.698681 
-2.888116 
-2.889404 
-3.060575 
-3.313477 
-3.420978 
-3.464977 
-3,758911 
-3.898113 
-3.965236 
-4.034424 
-4.373697 
-4.382095 
-4.390918 
-4.40P008 
-4.541711 
-4.632813 
-4.729350 
-4.737253 
-4.835818 
-6.129150 
-8.653698 
-8.787143 
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DISCRIMINANT 

ft  s 4.303303 

EVENT  StiQ  * 

EVENT  DESGN  ft 

D1SCR,  VALUE 

31 .00033 

.11.00000 

0.8952192 

7.003000 

18.00000 

0.7463280 

8.030000 

19.00000 

0.5422524 

35.00000 

53.00003 

0.5049603 

10.00000 

1,000000 

0.3952518 

9. 0i1?003 

20  » 30000 

0.464357 1E-0 1 

12.00000 

7.000000 

0.0000000 

29.00303 

17.03030 

-C.8439064E-01 

27.03000 

21.00000 

-0.9485781E-01 

24.00000 

149.0000 

-0.1727462 

2.030000 

46,00000 

-0.2409449 

30.00000 

16,00000 

-0.4555215 

26.00000 

77.00000 

-0.4670366 

11.30000 

3. 000000 

-0.5195667 

20.00000 

14,00000 

-0.6651423 

19.00000 

39,00000 

-0.8017970 

1.030000 

30,00000 

-0.8054416 

21.00000 

61,00000 

-0.8079039 

25.03000 

151.0000 

-0.B099695 

23.00000 

143.0000 

-0.9101455 

6.000000 

73.00000 

-0.9852136 

28.00000 

22.30000 

-l,05549«l 

5.O0PD00 

66.00000 

-1,072161 

17.03030 

49,00000 

-1.073619 

13.00000 

34.00000 

‘ -1.082644 

18.00000 

50.00000 

-1.113121 

3.000000 

55.00300 

-1.134726 

'4. 303000 

56.00020 

-1,324359 

33.0003W 

36.00000 

-1,345183 

16.00000 

39.00333 

-1,410806 

15,03030 

37.03000 

-1.549138 

22.03000 

62,00000 

-1.646932 

14.00030 

35.00000 

-1,710390 

34.00000 

47.00300 

-1.996634 

32.00000 

30.O0900 

-2.421961 
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SORTED 

DISCRIMINANT 

event  sed  tt 

10.00000 

19.0- 1000 
12. 30300 
24.  001400 

18.00000 

30.00000 
9,(303003 

34.00000 

32.000- 00 

25.00300 

22.00000 
14.0000(3 

I. 030000 
27.03000 
16.00000 

13.00000 
3.330300 

23.00300 

b.000000 

4.000330 

15.00000 

II. 00003 
3 1 .00013 
?.  1.000  JO 
35,0037(1 
33.0:1030 

7.000000 
23.OM0O0 
29.03300 
28.00030 
5.033033 
8.300300 

17.00300 

2.000000 


VENf  OISCRIMIN 
0=  5.000003 

EVENT  DESGN  « 

1.000003 

59.00000 

7.000000 
144.3330 
53.33300 

16.00030 

20.00000 

77.03000 
47.03003 
30,03303 

151.0003 

62.33000 

35.00000 

38.03000 

21.03000 

39.00000 

34.33000 

55.00000 

Id. 00000 

73.00000 

56.30000 

37.00000 

3.000003 

81.03000 
61.00300 
53.20000 
36. 0300? 

18.33000 

143.3000 

17.00000 

22.00000 

68.30000 

19.00003 

49.00030 
48,03033 


DlSCR.  VALUE 
0.5145466 
0.4549572 
0.3399252 
0.3394221 
0.2319742 
0,2078829 
0.1871995 
0. 1848801 
0.1664774 
0.1616617 
0.1409117 
0.1335987 
0.825C213E-01 
0.3482962E-01 
0.1236554E-01 
0.1154971E-01 
0.6297946E-02 
-0.1206172E-01 
-0.1532853E-01 
-0.2377582E-01 
-0.3080845E-01 
-0.3610981E-01 
-0.4606116E-01 
-0.5767143E-01 
-0.6189013E-01 
-0.6466675E-01 
-O.8e50597E-01 
*0.9410161 E-0 1 
-0.9704590E-01 
-0.1195118 
-0.1367149 
-0.1625949 
-0.1724290 
-0.2161083 
-0.3064649 


TABLE  V- 5 


EVENTS  SORTED 


SORTED 

DISCRIMINANT 

event  sea  u 

14.00030 

33.00000 

15.00000 
13.01000 
16.0P00W 

22.00000 
25.000KO 
<1.000000 

3.000000 

26.00000 
20.00000 
2 « 000300 

1.000030 
21.00000 

5.000000 

6.000000 

19.00000 

23.00000 
3H.03J00 

17.00000 
?B .03000 

32.00000 

11.00000 
ft. 000030 
34.3O000 
12.00000 

31.00000 

24.00000 

29.00000 

27.00000 
7.003000 

18.00000 
9.003300 
10.00000 
.35.30030 
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EVENT  DISCRIMINANTS 


As  6.000000 

EVENT  OESGN  * 

35.00000 

36.00000 

37.00000 

34.00000 
39,30000 

62.00000 

151.0000 

56.00000 
55.30G00 

77.00000 

14.00000 

40.00000 

30.00000 
61.0000e 

66.00000 

73.00000 

59.00000 

143.0000 

16.00000 

'49,00003 

22.00000 

.33,30000 

3.000000 

19.00000 

47.00000 
7.010000 

61.00000 

149.0300 

17.00300 

21,00000 

16.00000 

50.00000 

20.00000 

1 , 000000 

53.00000 


D1SCR.  VALUE 
1,199142 
1.194187 
1.163749 
1.019934 
0.9970392 
0.9834928 
0.9554608 
0.9206046 
0.8880939 
0.6239590 
0.7839437 
0,7839420 
0.7774760 
0.7498420 
0.6697786 
0.5812969 
0.5497030 
0.5256926 
0.3850524 
0.3744305 
8.3267765 
0.2169040 
0.1942062 
0.8777606E-01 
0.8527875E-01 
0.0000000 
-0.1198176 
-0.1295769 
-U. 1831255 
-0.2439444 
-0.2526371 
-0.3180960 
-0.3837500 
-3.4147766 
-1.033591 
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SORTED 

EVENT 

DISCRIMINANTS 

)I  SCR  J 1IN/VNT 

■»s 

7.0O0POO 

EVENT  seu  * 

EVENT  0E5GN  * 

D1SCR.  VALUE 

35.(10000 

53.0000? 

2.033586 

9.230000 

20.OOO00 

0.9826957 

10.03OO0 

1.000000 

0.7113269 

18.00000 

30. 33O0O 

0.6856962 

27. 00000 

21.00000 

0.8087300 

24,00300 

149,0300 

0.1970680 

29.T.1000 

17,00030 

0.  1868126 

31.OPO0O 

81. 30330 

-3.5148679E-01 

7.0O30OO 

18.30000 

-0.1062366 

28. 03000 

22,00000 

-0.1970615 

32.0O0W0 

30.00000 

-0.704H051 

0.000000 

19.00000 

-0.7690734 

30.00000 

47.00000 

-0.7986484 

19,03300 

59,00000 

-e. 899(653 

17,00000 

89,00000 

-1.115621 

6.003030 

T3. 00000 

-1.739066 

11.00000 

3.000030 

-1,843077 

23.00000 

193,0300 

-2.026750 

39. 00300 

16.00000 

-2,098837 

25.00000 

151.0000 

-2.188510 

23.00000 

14,00000 

-2.339602 

13,00000 

34.00000 

-2.368003 

26.03000 

77,00000 

-2.881839 

5.000030 

60,00000 

-2.530518 

14.00000 

35.00000 

-2.607703 

22.00303 

62.OOO0C? 

-2.705066 

16,03300 

39.OOOO0 

-2.796709 

21.00000 

61,00000 

-2.946625 

4.003003 

56.00000 

-2.967228 

1 .000000 

38.00000 

-3.353980 

3.000000 

55.00000 

-3.362595 

4. 3 00 OOP 

88.00000 

-3.368868 

33.0^300 

36,03000 

-3.473007 

15.02300 

37.30000 

-3.644625 

12.00003 

7.00000? 

-3.770361 
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SORTED  EVENT  DISCRIMINANTS 


DISCRIMINANT 
event  SEU  M 
10.00000 

7.000000 
27  • 00030 
31.DO0O0 
11.000(30 

24 .00000 

35.00000 

19.00000 

20.00000 

28.00000 
8.000000 

29.00000 

18.00000 

9.000000 

17.00000 

26.00000 

32.00000 

14.00000 

34.00000 

23.00000 

25.00000 

22.00000 

13.00000 

12.00000 

5.000000 

16.00000 
4.030030 
1.000000 
21.03000 
2.000030 

33.00000 

6.000000 

15.00000 
3.000300 

30.00000 


4:  8.000000 

EVENT  DESGN  ft 

1.000000 

18.00000 

21.00000 

81,00000 

3.000000 

149.0000 

53.00000 

59.00000 

14.00000 

22.00000 

19.00000 

17.00000 

50.00000 

20.00000 

49.00000 

77.00000 

30.00000 
35.(30000 

47.00000 

143.0000 

151.0000 

62.00000 

34.00000 
7.(300000 

68.00000 

39.00000 

56.00000 

38.00000 

61. 00000 

48.00000 

36.00000 

73.00000 

37.00000 

55.00000 
16.33800 


D1SCR.  VALUE 
3.3370O9 
2.062326 
-0.987O589C-01 
-0.3415936 
-0.4250169 
-0,6374735 
-0, 7380533 
-1,197331 
-1,416382 
-1,417003 
-1,459015 
-1,499571 
-1.532145 
-1.646327 
-1,814919 
-2,087240 
-2,188503 
-2.206860 
-2.354964 
-2.379538 
-2,470126 
-2,523605 
-2.537966 
-2,638682 
-2.989014 
-3.007554 
-3.102294 
-3.247847 
-3.457743 

-3.486869 

-3.787739 

-4.011271 

-4.221710 

-4.243390 

-4.306190 
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SORTED 

EVENT  DISCRIMINANTS 

DISCRIMINANT 

«s  9,000000 

KVKf-T  SEG  U 

EVENT  DESCN  9 

DISCR.  VALUE 

31.00300 

Rl, 30000 

-1.355317 

0.090300 

19.O0O00 

-1.909693 

7.900000 

18.00O0O 

-1.502222 

10.00900 

1,000000 

-1.565203 

35.00000 

53.30000 

-1.652115 

30.30000 

16.00000 

-1.781093 

26.00000 

77.00000 

-1,886929 

19.00000 

59.30000 

-2.139279 

29.00000 

199.0000 

-2.181797 

10.00000 

39.00O0O 

-2,393639 

12.3M300 

T. J0OO00 

-2,679137 

18.03003 

53.00000 

-2.723585 

29.00000 

17,00000 

-2.796755 

25.00000 

151.0000 

-2,839909 

23.00300 

193.0000 

-2,857876 

32.00000 

33.O0000 

-2,900297 

3.000080 

55,00000 

-2.915809 

21.03000 

61 .00000 

-2.960139 

9.000000 

56.00000 

-3,033518 

17.33000 

99.00000 

-3.083925 

6.000000 

73,00000 

-3.113988 

15.00000 

37.00000 

-3.165959 

20.00000 

19.33O00 

-3.297396 

1.000003 

38.00000 

-3.930501 

27.00000 

21.00000 

-3.505096 

9.000000 

20,00000 

-3.513059 

2.000000 

98.00000 

-3.513917 

33.00030 

36.00000 

-3,550718 

?fl. 0*030 

22,30000 

-3.568392 

5.003000 

68.00000 

-3.722703 

13.00000 

39,00000 

-9.121396 

22.31030 

o2, 09000 

-9.730292 

11.03300 

3.030000 

-9.851397 

19,00003 

35.00000 

-9,893920 

39.00000 

97,03000 

-5.610709 
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SORTED  EVENT  DISCRIMINANTS 


DlSr.RIMINANi  * 
EVEN!  SEJ  # 
20.P300O 
29.02000 
7.(300023 
30.03033 
8.035003 

35.00000 

20.00000 
12.33000 
31.00300 

27 . 00000 
22.03003 

9. 030000 

1 9.00000 

29.03000 

11.00000 
(> . 003000 
10.0(3000 

16.03000 
13.003110 

23.03000 

33.00000 

1.000000 

9.000000 

5. 000000 

2.000000 

3.000000 

21.00  VO 

15.03000 
>4.03000 
14.00030 

10.30000 

32.00000 

26.00000 

1 7, 00.130 

25.0O0OO 


10.00000 
1VENT  PESGN  * 

14.00000 
17.001*'  ‘ 

16.00000 
16.00000 
19.00030 

53.00000 
22.O000P 

7.000000 

61.00000 
21.00000 

62.30000 
20.0UOO0 
59,30030 

149.0000 

3.000000 
73.P0O0O 

50.30000 

39.00000 

34.00000 

143.0000 

36.00000 

30.00000 

56.00000 

60.00000 

40.00000 

55.00000 

61.30000 

37.00000 

47.00000 
35.00003 

1.000000 

30.00000 

77.00000 

49.00000 

151.0000 


D1SCR.  VALUE 
0.8262868 
0.0843088 
0.5886557 
0.5604186 
0.5275884 
0.4809289 
0.3926786 
0,3803024 
0,3647542 
0.3487525 
0,2880888 
0.2545587 
0.2498338 
0.6007352E-01 
0.5403280E-01 
0.4966974E-01 
0.4880795E-01 
U.4675579E-01 
0. 1254082E-01 
0.5267709E-02 
0.0000000 
0.0000000 
0.0000000 
0.0000000 
0.0000000 
0.0000000 
0,0000000 
0.0000000 
-0.5356172E-31 
-0.1132622 
-0.1137602 
-0.1579071 
-0.1667642 
-0,2182080 
-e. 2321825 
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EVENTS  SORTED  BY  MEASURED  DISCRIMINANT  VALUE 

(PAGE  11  OF  13) 


SORTED  CVCNT  DISCRIMINANTS 
DISCRIMINANT  Js  11.00000 
EVENT  SHU  w EVENT  DESGN  » DISCR.  VALUE 
2O.O0OOU  22.0O3O0  2.2930/13 

23.00000  19.00000  2,105371 

6.OOO0PO  73.00000  2.072972 

10.30900  1,009000  1.90859a 

13.00000  39.09000  1.919897 

22.09000  62. 00000  1. 919817 

11.00000  3.000000  1.891651 

25.0000P  151.0000  1.859027 

19.00000  35.00000  1.755218 

7.000000  18,03000  1,679335 

27.00090  21.00000  1.652901 

29.030HO  17.00D00  1.637867 

9.000000  23.00300  1.583375 

19.00090  59.O0O30  1.553226 

26.30030  77.00000  1.532199 

8.030000  19.00000  1.502919 

12.00000  7 , 000O00  1.976109 

35.00000  53,00000  1,972136 

30.00090  16.00000  1.956291 

23.00000  193.0000  1.913250 

16.03000  39.00000  1.395308 

18.00000  50.0OOO0  1.325080 

32.00000  3O.00000  1,280355 

29.00000  199.0000  1,279927 

31 .00000  81.00000  1.251103 

17.00000  99.00000  1.219591 

39 ,00030  97.00000  1.195368 

33.00000  36,00000  0.0000000 

1.000000  38,03000  0,0000000 

9.000000  56.00030  0.0000000 

2.000000  98. 110000  0.0&000PO 

3.000000  55.00000  0.0000000 

21.33000  61.00000  0.0000000 

15.00000  37.00000  0.0000000 

5.0009100  68.00000  0,0000000 
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EVENTS  SORTED  BY  MEASURED  DISCRIMINANT  VALUE 

(PAGE  12  OF  13) 


SORTED  EVF'JT  DISCRIMINANTS 
niSC.RlMINANT  «s  12.00000 
EVENT  sea  » EVE. a desgn  # piscr.  VALUE 


14.00000 

12.00300 
16.OU3O0 

11.03000 
6.0000UU 
19.33000 

20.30000 

18.00000 

26.00000 
2 '<.00000 
8.00000k) 

28.00000 

29.00000 

35.03000 

10.00000 

23.00000 

34.03000 

7.000000 

25.00000 

27.00000 

31.00000 

9.000000 

17.00000 
30. JOJ30 

32.00000 

33.0O03O 

1.000000 
'4. 0OP000 

5.000000 

2.000000 
3.30OOO0 
22.00003 
21.30300 

13.30000 

15.30000 


35.00000 
7.U0000H 

39.30000 
3.033000 

73.00300 

59.00000 

16.00000 

50 .00000 

77.00000 

169.0000 
19.00PO0 

22.30000 
17.P0O00 

53.00000 

1.000000 

163.0300 

67.03000 

18.00030 

151.0000 

21.00000 

81.00000 
2R.03P0O 
69.09000 

16.00030 

30.00000 

36.03000 

38.00300 

56.00000 
68.90000 

66.00000 

55.00000 

62.00000 

61.00000 

36.00000 

37.00000 


2.380652 
1,033109 
1.768305 
1.760996 
1.733068 
1.590322 
1.556960 
1.551296 
1.659213 
1.625363 
1,396592 
1.391832 
1.389096 
1.358329 
1,346107 
1,322733 
1.316248 
1.313726 
1,305879 
1.305663 
1.294377 
1.279096 
1.221704 
1.146245 
1.064979 
0.0000000 
0.0000000 
0.0000000 
0.0900000 
0.0000000 
0 . 0000000 
0.0000000 
0.0000000 
0.0000000 
0,0000000 
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EVENTS  SORTED  BY  MEASURED  DISCRIMINANT  VALUE 

(PAGE  13  OF  13) 


SORTED  EVENT  DISCRIMINANTS 


MULTIVARIATE 

DISCRIMINANT 

EVENT  DESGN  # 

DISCR,  VALUE 

ie. 00000 

91,92015 

53.003HC 

72.91619 

17.00000 

71.10209 

19.00000 

68.32903 

1 9.0:10  J9 

66.60777 

81.00000 

62.26123 

lo.POUU'? 

56.76215 

21.00000 

51.16705 

1.000000 

96.33960 

22.00000 

38,05307 

20.00000 

37.93552 

59.00C0P 

39,  12929 

7.00OO00 

30,59122 

1 -T9, 0000 

28.20879 

o2, 00000 

23.57536 

3.000000 

23.57518 

50.00000 

9.972989 

77.00000 

9.228659 

'*8.00000 

6.897387 

97.00000 

3.655029 

39.00000 

3.212139 

99.00000 

2.599737 

30.00000 

1.909399 

36.00000 

1.501603 

38.00000 

1.993665 

ot ,00000 

1.220119 

73, 00O JO 

0. 1)000000 

66,00000 

0. 0000000 

55,0000? 

0.0000000 

56, 00300 

3.0300000 

151.0000 

0.0000000 

37.00000 

0.0000000 

193,0000 

0.0000000 

39,00000 

0.0000000 

35,00000 

0.0000000 
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TABLE  V-  6 

EVENT  IDENTIFICATION  RESULTS 
(PACE  1 OF  2) 


INITIAL  NORMAL  EARTHQUAKE  TRAINING  SET 

Event 

1 

Event 

Comments 

Sequence  Number 

Designation  Number 

5 

68 

Dropped  due  to  program  error. 

6 

73 

11 

3 

13 

34 

14 

3.5 

16 

35 

23 

143 

24 

149 

25 

151 

26 

77 

32 

30 

34 

47 

UNKNOWN  EVENTS  CLASSIFIED  AS  EARTHQUAKES 

Event 

Event 

Comments 

Sequence  Number 

Designation  Number 

3 55* 

4 56 

15  37 

21  61 

33  36 

17  49 
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EVENT  IDENTIFICATION  RESULTS 
(PAGE  2 OF  2) 


CORRELATED  EVENT  GROUP  NUMBER  ONE 

Event 

Event 

Comments 

Sequence  Number 

Designation  Number 

7 

18 

8 

19 

20 

14 

29 

17 

9 

20 

10 

1 

12 

7 

Dropped  due  to  program  error. 

19 

59' 

22 

62 

27 

21 

28 

22 

30 

16 

31 

81 

35 

53 

UNKNOWN  EVENTS  NEITHER  CLASSIFIED  AS  GROUP  NUMBER  ONE  OR 

AS  NORMAL  SHALLOW  FOCUS  EARTHQUAKES 

Event 

1 

Event 

Comments 

Sequence  Number 

Designation  Number 

-25 


TABLE  V-7 
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DETECTABILITY  OF  CORRELATED  EVENT  GROUP  NUMBER  ONE 
AND  UNCORRELATED  UNKNOWNS 


Discriminant 

Number 

Detectability 

Group  No.  One 

Event  38 

Event  48 

Event  50 

1 

0.  94 

5.  05 

1. 05 

-0.  90 

2 

1.13 

1.  32 

1. 86 

-1.  34 

3 

0.  96 

1. 46 

0.  71 

-2.  50 

4 

2.  37 

0.  65 

1.  99 

-0.  07 

5 

-0.  05 

0.  75 

-4.  93 

4.  00 

6 

-1. 54 

2.  62 

2.  62 

0.  82 

7 

1. 61 

-0.  71 

-0.  71 

2.  15 

8 

2.  02 

-1.  50 

-1.  96 

1.  81 

9 

2. 12 

-0.  95 

-1.  14 

0.  74 

10 

5.  59 

N.  A.* 

N.  A. 

1. 00 

11 

0.  40 

N.  A. 

N.  A. 

-0.  35 

12 

-0.41 

N.  A. 

N.  A. 

0.  71 

MVD 

9.  62 

-0.  95 

1. 07 

1. 67 

N.  A.  : Not  available  due  to  lack  of  detections. 


identification  process  will  retain  these  events  as  unknown  unusual  events  un- 
til sufficient  event  data  are  entered  into  the  system  to  group  them  into  an- 
other correlated  grouping  of  similar  events. 

For  robustness,  the  detectabilities  were  estimated  by  comput- 
ing the  difference  between  the  abnormal  group  or  abnormal  event  median  and 
the  median  of  the  updated  earthquake  population.  Abnormal  event  medians 
are  taken  as  the  observed  discriminant  value.  This  difference  of  medians 
was  divided  by  the  normal  earthquake  population  deviation;  computed  as  the 
difference  between  the  85%  level  of  normal  earthquakes  and  the  median  of 
normal  earthquakes  as  shown  on  Figure  V-l. 


The  detectability  of  correlated  group  number  one  for  discrimi- 
nants number  one  through  twelve  were  used  as  weights  for  combining  discrimi- 
nant observations  into  a single  multivariate  discriminant.  The  operator  for 
applying  these  group  projection  weights  to  an  observed  set  of  discriminants 
is  shown  by  equation  (IV-1).  The  results  of  the  multivariate  discriminant 
calculations  sorted  by  discriminant  value  are  shown  on  part  13  of  Table  V-5. 
Medians  and  deviations  are  shown  on  Tables  V-2  through  V-4.  Remember 
that  each  observation  of  an  unknown  event  discriminant  vector  is  transformed 
into  a detectability  vector  by  subtracting  the  normal  earthquake  median  and 
dividing  by  the  normal  earthquake  deviation.  This  converts  each  discrimi- 
nant observation  vector  into  a homogeneous  set  of  dimensionless  components 
of  a z-vector.  To  this,  the  projection  operator  is  applied  with  respect  to  the 
correlated  group  z-vector.  A projection  above  a prescribed  threshold  is  ac- 
cepted as  a new  similar  member  of  the  group. 
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SECTION  VI 
CONCLUSIONS 


Our  goal  was  to  develop  a package  to  objectively  evaluate  the 
performance  of  a set  of  selected  discriminants  as  potential  components  of  a 
seismic  event  identification  system.  As  a ground  rule,  the  following  opera- 
tions were  performed  automatically  to  achieve  an  objective  evaluation  of  per- 
formance. 

• Determine  ground  motion  by  cr  rrecting  for  system  response. 

• Edit  phases. 

• Extract  waveforms  and  determine  detection  status. 

• Filter  waveforms. 

• Generate  a signal  measurement  file. 

• Generate  an  event  source  measurement  file  by  operating  on 

the  signal  measurement  file. 

• Generate  an  event  discriminant  file  by  operating  on  the  source 
measurement  file. 

• Generate  inductive  grouping  of  event  types  by  operating  on  the 
discriminant  file. 

The  final  results  of  this  mode  of  discriminant  processing  are  reproducibly 
determined  subject  to  specification  of  a set  of  algorithms  for  determining  the 
detection  status  of  edits,  algorithms  for  measuring  signals  or  noise  in  edit 
windows,  algorithms  for  converting  signal  measurements  to  unbiased  esti- 
mates of  source  parameters,  algorithms  for  converting  deriving  discriminants 
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from  source  parameter  measurements,  and  decision  functions  for  grouping 
events  of  similar  discriminant  characteristics.  In  addition,  to  start  the  pro- 
cess it  was  necessary  to  introduce  a training  set  of  normal  shallow  earth- 
quakes, determine  the  statistical  characteristics  of  their  discriminants  and 
use  those  characteristics  to  normalize  the  observed  discriminants  of  subse- 
quent unknown  events. 

Our  initial  test  run  was  performed  on  thirty  five  events,  twelve 
of  which  were  selected  as  the  initial  normal  shallow  earthquake  training  set. 

Of  the  unknown  set,  six  were  correlated  with  normal  shallow  earthquakes  and 
fourteen  were  correlated  with  a similar  group  of  events  which  are  not  normal 
shallow  earthquakes,  and  three  remained  as  unknown  unusual  events  which 
were  neither  classified  as  normal  shallow  earthquakes  nor  as  the  correlated 
group  of  unusual  events.  This  test  was  performed  with  arbitrarily  selected 
decision  thresholds  for  grouping  events.  Based  on  statistics  determined 
from  the  test  and  the  observed  grouping  of  events  a threshold  can  be  set  to 
operate  with  any  desired  operating  characteristics.  For  example,  two  normal 
earthquakes  appeared  to  be  included  in  the  correlated  group.  By  raising  the 
threshold,  the  two  apparent  normal  earthquakes  could  have  been  classified 
correctly  without  removing  any  other  events  from  the  correlated  group.  Two 
of  the  unknown  unusual  events  were  deep  events  from  the  Afghan- USSR  border 
which  appear  to  be  similar  to  each  other.  The  third  unusual  event  was  an 
intermediate  depth  event  from  the  Kurile  Islands. 

Relative  ranking  performance  of  particular  discriminants  are 
apparent  by  examination  of  tables  of  likelihood  tests  and  detectabilities  shown 
in  the  preceding  section.  It  should  be  kept  in  mind,  however,  that  such  rank- 
ings of  discriminant  effectiveness  are  only  relevant  to  the  single  correlated 
gr  :up  thus  far  obtained.  Other  apparently  less  effective  discriminants  may 
be  important  in  detecting  members  of  other  unusual  event  groupings  not  yet 
observed.  Once  the  decision  rules  are  set,  the  package  provides  a multi-user 


environment.  The  classification  files  can  be  used  as  permanent,  maintainable 
and  transportable  files.  As  more  events  are  added  into  such  a system,  the 
performance  will  improve  as  a result  of  building  up  the  grouped  populations  of 
abnormal  events  and  normal  earthquake  populations.  This  will  result  in  better 
statistical  estimates  of  population  characteristics.  In  another  sense,  better 
performance  can  be  anticipated  by  expanding  the  list  of  discriminants  and  the 
number  of  stations  in  the  network. 
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APPENDIX  A 


VFM  NARROWBAND  FILTER  DESCRIPTION 


A cusp- shaped  narrowband  filter  is  used  in  many  different  ap- 
plications throughout  the  routines  used  in  obtaining  discriminant  measure- 
ments. This  filter  gives  optimal  time  and  frequency  domain  characteristics 
within  the  limits  of  the  Uncertainty  Principle.  The  goals  which  this  filter 
attempts  to  satisfy  are:  (1)  minimum  width  in  the  frequency  domain,  and  (2) 
maximum  ripple  suppression  in  the  time  domain.  Since  the  Uncertainty  Prin- 
ciple disallows  the  possibility  of  simultaneously  satisfying  these  two  goals  to 
an  arbitrary  precision,  the  filter  was  chosen  to  optimize  these  characteris- 
tics. This  filter  has  been  used  to  obtain  accurate  time  separation  and  ac- 
curate amplitudes  in  separating  multiple  explosions.  It  has  also  been  used 
in  calculating  variable  frequency  magnitudes  (VFM)  based  on  the  maximum 
envelope  amplitudes  for  the  filtered  signal. 

Figure  A-l  illustrates  the  filter  function  in  the  frequency  do- 
main. The  cusp- shaped  filter,  F,  is  of  the  form: 


where 


ff  - (f  - 3/2  Af) I 

1 - cos  tt/3  I — , (f  - 3/2Af)  < f < f 

i Af  J c c 

r (f  + 3/2  Af)  - fl 

1 - cos  TT/3  — — , f < f < (f  + 3/2 Af) 

L Af  J C C 

0 , otherwise 

f = frequency 

f = frequency  at  which  the  filter  is  centered 

Af  = the  desired  width  (in  Hz)  at  \ amplitude  points  (-6dB  points). 


FIGURE  A- 1 

NARROWBAND  FILTER  USED  FOR  VFM  MEASUREMENT 
THE  WIDTH  AT  ONE-HALF  MAXIMUM  AMPLITUDE 
IS  DESIGNATED  Af 
(Lambert,  Bache,  and  Savino,  1977) 
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